source: trunk/exports/G2export_CIF.py @ 5127

Last change on this file since 5127 was 5127, checked in by toby, 9 months ago

misc CIF fixes (init block name; stop save of default template file name)

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 248.0 KB
Line 
1#!/usr/bin/env python
2# -*- coding: utf-8 -*-
3########### SVN repository information ###################
4# $Date: 2022-01-07 21:17:14 +0000 (Fri, 07 Jan 2022) $
5# $Author: toby $
6# $Revision: 5127 $
7# $URL: trunk/exports/G2export_CIF.py $
8# $Id: G2export_CIF.py 5127 2022-01-07 21:17:14Z toby $
9########### SVN repository information ###################
10'''
11*Module G2export_CIF: CIF Exports*
12------------------------------------------------------
13
14This implements a complex exporter :class:`ExportCIF` that can implement an
15entire project in a complete CIF intended for submission as a
16publication. In addition, there are three subclasses of :class:`ExportCIF`:
17:class:`ExportProjectCIF`,
18:class:`ExportPhaseCIF` and :class:`ExportDataCIF` where extra parameters
19for the _Exporter() determine if a project, single phase or data set are written.
20'''
21
22from __future__ import division, print_function
23import platform
24import datetime as dt
25import os.path
26import sys
27import numpy as np
28if '2' in platform.python_version_tuple()[0]:
29    import cPickle as pickle
30else:
31    import pickle
32import copy
33import re
34try:
35    import wx
36    import wx.lib.scrolledpanel as wxscroll
37    import wx.lib.resizewidget as rw
38except ImportError:
39    # Avoid wx dependency for CLI
40    class Placeholder(object):
41        def __init__(self):
42            self.BoxSizer = object
43            self.Button = object
44            self.Dialog = object
45            self.ScrolledPanel = object
46    wx = Placeholder()
47    wxscroll = Placeholder()
48import GSASIIpath
49GSASIIpath.SetVersionNumber("$Revision: 5127 $")
50import GSASIIIO as G2IO
51try:
52    import GSASIIctrlGUI as G2G
53except ImportError:
54    pass
55import GSASIIobj as G2obj
56import GSASIImath as G2mth
57import GSASIIspc as G2spc
58import GSASIIlattice as G2lat
59import GSASIIstrMain as G2stMn
60import GSASIIstrIO as G2stIO       
61import GSASIImapvars as G2mv
62import GSASIIElem as G2el
63import GSASIIpy3 as G2py3
64
65DEBUG = False    #True to skip printing of reflection/powder profile lists
66
67CIFdic = None
68
69cellNames = ['length_a','length_b','length_c',
70             'angle_alpha','angle_beta ','angle_gamma',
71             'volume']
72def striphist(var,insChar=''):
73    'strip a histogram number from a var name'
74    sv = var.split(':')
75    if len(sv) <= 1: return var
76    if sv[1]:
77        sv[1] = insChar
78    return ':'.join(sv)
79
80def getCellwStrain(phasedict,seqData,pId,histname):
81    'Get cell parameters and their errors for a sequential fit'
82    #newCellDict = {}
83    #if name in seqData and 'newCellDict' in seqData[histname]:
84    #    newCellDict.update(seqData[histname]['newCellDict'])
85   
86    pfx = str(pId)+'::' # prefix for A values from phase
87    Albls = [pfx+'A'+str(i) for i in range(6)]
88    Avals = G2lat.cell2A(phasedict['General']['Cell'][1:7])
89    #AiLookup = {}
90    DijLookup = {}
91    zeroDict = dict(zip(Avals,6*[0.,]))
92    for i,v in enumerate(('D11','D22','D33','D12','D13','D23')):
93        if pfx+v in seqData[histname]['newCellDict']:
94            Avals[i] = seqData[histname]['newCellDict'][pfx+v][1]
95            #AiLookup[seqData[histname]['newCellDict'][pfx+v][0]] = pfx+v
96            DijLookup[pfx+v] = seqData[histname]['newCellDict'][pfx+v][0]
97    covData = {  # relabeled with p:h:Dij as p::Ai
98        'varyList': [DijLookup.get(striphist(v),v) for v in seqData[histname]['varyList']], 
99        'covMatrix': seqData[histname]['covMatrix']}
100    # apply symmetry
101    cellDict = dict(zip(Albls,Avals))
102    try:    # convert to direct cell
103        A,zeros = G2stIO.cellFill(pfx,phasedict['General']['SGData'],cellDict,zeroDict)
104        cell = list(G2lat.A2cell(A)) + [G2lat.calc_V(A)]
105        cE = G2stIO.getCellEsd(pfx,phasedict['General']['SGData'],A,covData)
106    except:
107        cell = 7*[None]
108        cE = 7*[None]
109    return cell,cE
110
111def mkSeqResTable(mode,seqHistList,seqData,Phases,Histograms,Controls):
112    '''Setup sequential results table (based on code from
113    GSASIIseqGUI.UpdateSeqResults)
114
115    TODO: This should be merged with the table build code in
116    GSASIIseqGUI.UpdateSeqResults and moved to somewhere non-GUI
117    like GSASIIstrIO to create a single routine that can be used
118    in both places, but this means returning some
119    of the code that has been removed from there
120    '''
121
122    newAtomDict = seqData[seqHistList[0]].get('newAtomDict',{}) # dict with atom positions; relative & absolute
123    atomLookup = {newAtomDict[item][0]:item for item in newAtomDict if item in seqData['varyList']}
124    phaseLookup = {Phases[phase]['pId']:phase for phase in Phases}
125
126    # make dict of varied cell parameters equivalents
127    ESDlookup = {} # provides the Dij term for each Ak term (where terms are refined)
128    Dlookup = {} # provides the Ak term for each Dij term (where terms are refined)
129    newCellDict = {}
130    for name in seqHistList:
131        if name in seqData and 'newCellDict' in seqData[name]:
132            newCellDict.update(seqData[name]['newCellDict'])
133    cellAlist = []
134    for item in newCellDict:
135        cellAlist.append(newCellDict[item][0])
136        if item in seqData.get('varyList',[]):
137            ESDlookup[newCellDict[item][0]] = item
138            Dlookup[item] = newCellDict[item][0]
139    # add coordinate equivalents to lookup table
140    for parm in atomLookup:
141        Dlookup[atomLookup[parm]] = parm
142        ESDlookup[parm] = atomLookup[parm]
143
144    # get unit cell & symmetry for all phases & initial stuff for later use
145    RecpCellTerms = {}
146    SGdata = {}
147    uniqCellIndx = {}
148    #initialCell = {}
149    RcellLbls = {}
150    zeroDict = {}
151    for phase in Phases:
152        pId = Phases[phase]['pId']
153        pfx = str(pId)+'::' # prefix for A values from phase
154        RcellLbls[pId] = [pfx+'A'+str(i) for i in range(6)]
155        RecpCellTerms[pId] = G2lat.cell2A(Phases[phase]['General']['Cell'][1:7])
156        zeroDict[pId] = dict(zip(RcellLbls[pId],6*[0.,]))
157        SGdata[pId] = Phases[phase]['General']['SGData']
158        laue = SGdata[pId]['SGLaue']
159        if laue == '2/m':
160            laue += SGdata[pId]['SGUniq']
161        for symlist,celllist in G2lat.UniqueCellByLaue:
162            if laue in symlist:
163                uniqCellIndx[pId] = celllist
164                break
165        else: # should not happen
166            uniqCellIndx[pId] = list(range(6))
167           
168    # scan for locations where the variables change
169    VaryListChanges = [] # histograms where there is a change
170    combinedVaryList = []
171    firstValueDict = {}
172    vallookup = {}
173    posdict = {}
174    prevVaryList = []
175    foundNames = []
176    missing = 0
177    for i,name in enumerate(seqHistList):
178        if name not in seqData:
179            if missing < 5:
180                print(" Warning: "+name+" not found")
181            elif missing == 5:
182                print (' Warning: more are missing')
183            missing += 1
184            continue
185        foundNames.append(name)
186        maxPWL = 5
187        for var,val,sig in zip(seqData[name]['varyList'],seqData[name]['variables'],seqData[name]['sig']):
188            svar = striphist(var,'*') # wild-carded
189            if 'PWL' in svar:
190                if int(svar.split(':')[-1]) > maxPWL:
191                    continue
192            if svar not in combinedVaryList:
193                # add variables to list as they appear
194                combinedVaryList.append(svar)
195                firstValueDict[svar] = (val,sig)
196        if prevVaryList != seqData[name]['varyList']: # this refinement has a different refinement list from previous
197            prevVaryList = seqData[name]['varyList']
198            vallookup[name] = dict(zip(seqData[name]['varyList'],seqData[name]['variables']))
199            posdict[name] = {}
200            for var in seqData[name]['varyList']:
201                svar = striphist(var,'*')
202                if 'PWL' in svar:
203                    if int(svar.split(':')[-1]) > maxPWL:
204                        continue
205                posdict[name][combinedVaryList.index(svar)] = svar
206            VaryListChanges.append(name)
207    if missing:
208        print (' Warning: Total of %d data sets missing from sequential results'%(missing))
209
210    #### --- start building table
211    histNames = foundNames
212#            sampleParms = GetSampleParms()
213    nRows = len(histNames)
214    tblValues = [list(range(nRows))]   # table of values arranged by columns
215    tblSigs =   [None]                 # a list of sigma values, or None if not defined
216    tblLabels = ['Number']             # a label for the column
217    tblTypes = ['int']
218    # start with Rwp values
219    tblValues += [[seqData[name]['Rvals']['Rwp'] for name in histNames]]
220    tblSigs += [None]
221    tblLabels += ['Rwp']
222    tblTypes += ['10,3']
223    # add changing sample parameters to table
224    sampleParmDict = {'Temperature':[],'Pressure':[],'Time':[],
225                 'FreePrm1':[],'FreePrm2':[],'FreePrm3':[],'Omega':[],
226                 'Chi':[],'Phi':[],'Azimuth':[],}
227    for key in sampleParmDict:
228        for h in histNames:
229            var = ":" + str(Histograms[h]['hId']) + ":" + key
230            sampleParmDict[key].append(seqData[h]['parmDict'].get(var))
231        if not np.all(np.array(sampleParmDict[key]) == sampleParmDict[key][0]):
232            tblValues += [sampleParmDict[key]]
233            tblSigs.append(None)
234            if 'FreePrm' in key and key in Controls:
235                tblLabels.append(Controls[item])
236            else:
237                tblLabels.append(key)
238            tblTypes += ['float']
239
240    # add unique cell parameters 
241    if len(newCellDict):
242        for pId in sorted(RecpCellTerms):
243            pfx = str(pId)+'::' # prefix for A values from phase
244            cells = []
245            cellESDs = []
246            Albls = [pfx+'A'+str(i) for i in range(6)]
247            for name in histNames:
248                #if name not in Histograms: continue
249                hId = Histograms[name]['hId']
250                phfx = '%d:%d:'%(pId,hId)
251                esdLookUp = {}
252                dLookup = {}
253                for item in seqData[name]['newCellDict']:
254                    if phfx+item.split('::')[1] in seqData[name]['varyList']:
255                        esdLookUp[newCellDict[item][0]] = item
256                        dLookup[item] = newCellDict[item][0]
257                covData = {'varyList': [dLookup.get(striphist(v),v) for v in seqData[name]['varyList']],
258                    'covMatrix': seqData[name]['covMatrix']}
259                A = RecpCellTerms[pId][:] # make copy of starting A values
260                # update with refined values
261                for i,j in enumerate(('D11','D22','D33','D12','D13','D23')):
262                    var = str(pId)+'::A'+str(i)
263                    Dvar = str(pId)+':'+str(hId)+':'+j
264                    # apply Dij value if non-zero
265                    if Dvar in seqData[name]['parmDict']:
266                        parmDict = seqData[name]['parmDict']
267                        if parmDict[Dvar] != 0.0:
268                            A[i] += parmDict[Dvar]
269                    # override with fit result if is Dij varied
270                    if var in cellAlist:
271                        try:
272                            A[i] = seqData[name]['newCellDict'][esdLookUp[var]][1] # get refined value
273                        except KeyError:
274                            pass
275                # apply symmetry
276                cellDict = dict(zip(Albls,A))
277                try:    # convert to direct cell
278                    A,zeros = G2stIO.cellFill(pfx,SGdata[pId],cellDict,zeroDict[pId])
279                    c = G2lat.A2cell(A)
280                    vol = G2lat.calc_V(A)
281                    cE = G2stIO.getCellEsd(pfx,SGdata[pId],A,covData)
282                except:
283                    c = 6*[None]
284                    cE = 6*[None]
285                    vol = None
286                # add only unique values to table
287                if name in Phases[phaseLookup[pId]]['Histograms']:
288                    cells += [[c[i] for i in uniqCellIndx[pId]]+[vol]]
289                    cellESDs += [[cE[i] for i in uniqCellIndx[pId]]+[cE[-1]]]
290                else:
291                    cells += [[None for i in uniqCellIndx[pId]]+[None]]
292                    cellESDs += [[None for i in uniqCellIndx[pId]]+[None]]
293            p = phaseLookup[pId]
294            tblLabels += ['{}, {}'.format(G2lat.cellAlbl[i],p) for i in uniqCellIndx[pId]]
295            tblTypes += ['10,5' if i <3 else '10,3' for i in uniqCellIndx[pId]]
296            tblLabels.append('{}, {}'.format('Volume',p))
297            tblTypes += ['10,3']
298            tblValues += zip(*cells)
299            tblSigs += zip(*cellESDs)
300
301    # sort out the variables in their selected order
302    varcols = 0
303    varlbls = []
304    for d in posdict.values():
305        varcols = max(varcols,max(d.keys())+1)
306    # get labels for each column
307    for i in range(varcols):
308        lbl = ''
309        for h in VaryListChanges:
310            if posdict[h].get(i):
311                if posdict[h].get(i) in lbl: continue
312                if lbl != "": lbl += '/'
313                lbl += posdict[h].get(i)
314        varlbls.append(lbl)
315    vals = []
316    esds = []
317    varsellist = None        # will be a list of variable names in the order they are selected to appear
318    # tabulate values for each hist, leaving None for blank columns
319    for name in histNames:
320        if name in posdict:
321            varsellist = [posdict[name].get(i) for i in range(varcols)]
322            # translate variable names to how they will be used in the headings
323            vs = [striphist(v,'*') for v in seqData[name]['varyList']]
324            # determine the index for each column (or None) in the seqData[]['variables'] and ['sig'] lists
325            sellist = [vs.index(v) if v is not None else None for v in varsellist]
326            #sellist = [i if striphist(v,'*') in varsellist else None for i,v in enumerate(seqData[name]['varyList'])]
327        if not varsellist: raise Exception()
328        vals.append([seqData[name]['variables'][s] if s is not None else None for s in sellist])
329        esds.append([seqData[name]['sig'][s] if s is not None else None for s in sellist])
330    tblValues += zip(*vals)
331    tblSigs += zip(*esds)
332    tblLabels += varlbls
333    tblTypes += ['float' for i in varlbls]
334   
335    # tabulate constrained variables, removing histogram numbers if needed
336    # from parameter label
337    depValDict = {}
338    depSigDict = {}
339    for name in histNames:
340        for var in seqData[name].get('depParmDict',{}):
341            val,sig = seqData[name]['depParmDict'][var]
342            svar = striphist(var,'*')
343            if svar not in depValDict:
344               depValDict[svar] = [val]
345               depSigDict[svar] = [sig]
346            else:
347               depValDict[svar].append(val)
348               depSigDict[svar].append(sig)
349
350    # add the dependent constrained variables to the table
351    for var in sorted(depValDict):
352        if len(depValDict[var]) != len(histNames): continue
353        tblLabels.append(var)
354        tblTypes.append('10,5')
355        tblSigs += [depSigDict[var]]
356        tblValues += [depValDict[var]]
357
358    # add refined atom parameters to table
359    for parm in sorted(atomLookup):
360        tblLabels.append(parm)
361        tblTypes.append('10,5')
362        tblValues += [[seqData[name]['newAtomDict'][atomLookup[parm]][1] for name in histNames]]
363        if atomLookup[parm] in seqData[histNames[0]]['varyList']:
364            col = seqData[histNames[0]]['varyList'].index(atomLookup[parm])
365            tblSigs += [[seqData[name]['sig'][col] for name in histNames]]
366        else:
367            tblSigs += [None]
368
369    # compute and add weight fractions to table if varied
370    for phase in Phases:
371        var = str(Phases[phase]['pId'])+':*:Scale'
372        if var not in combinedVaryList+list(depValDict.keys()): continue
373        wtFrList = []
374        sigwtFrList = []
375        for i,name in enumerate(histNames):
376            if name not in Phases[phase]['Histograms']:
377                wtFrList.append(None)
378                sigwtFrList.append(0.0)
379                continue
380            elif not Phases[phase]['Histograms'][name]['Use']:
381                wtFrList.append(None)
382                sigwtFrList.append(0.0)
383                continue
384            wtFrSum = 0.
385            for phase1 in Phases:
386                if name not in Phases[phase1]['Histograms']: continue
387                if not Phases[phase1]['Histograms'][name]['Use']: continue
388                wtFrSum += Phases[phase1]['Histograms'][name]['Scale'][0]*Phases[phase1]['General']['Mass']
389            var = str(Phases[phase]['pId'])+':'+str(i)+':Scale'
390            wtFr = Phases[phase]['Histograms'][name]['Scale'][0]*Phases[phase]['General']['Mass']/wtFrSum
391            wtFrList.append(wtFr)
392            if var in seqData[name]['varyList']:
393                sig = seqData[name]['sig'][seqData[name]['varyList'].index(var)]*wtFr/Phases[phase]['Histograms'][name]['Scale'][0]
394            elif var in seqData[name].get('depParmDict',{}):
395                _,sig = seqData[name]['depParmDict'][var]
396            else:
397                sig = 0.0
398            sigwtFrList.append(sig)
399        p = phaseLookup[Phases[phase]['pId']]
400        tblLabels.append(p + ' Wgt Frac')
401        tblTypes.append('10,4')
402        tblValues += [wtFrList]
403        tblSigs += [sigwtFrList]
404    return tblLabels,tblValues,tblSigs,tblTypes
405
406
407# Refactored over here to allow access by GSASIIscriptable.py
408def WriteCIFitem(fp, name, value=''):
409    '''Helper function for writing CIF output.'''
410    # Ignores unicode issues
411    if value:
412        if "\n" in value or (len(value) > 70 and ' ' in value.strip()):
413            if name.strip():
414                fp.write(name+'\n')
415            fp.write(';\n'+value+'\n')
416            fp.write(';'+'\n')
417        elif " " in value:
418            if len(name)+len(value) > 65:
419                fp.write(name + '\n   ' + '"' + str(value) + '"'+'\n')
420            else:
421                fp.write(name + '  ' + '"' + str(value) + '"'+'\n')
422        else:
423            if len(name)+len(value) > 65:
424                fp.write(name+'\n   ' + value+'\n')
425            else:
426                fp.write(name+'  ' + value+'\n')
427    else:
428        fp.write(name+'\n')
429
430def RBheader(fp):
431    WriteCIFitem(fp,'\n# RIGID BODY DETAILS')
432    WriteCIFitem(fp,'loop_\n    _restr_rigid_body_class.class_id\n    _restr_rigid_body_class.details')
433
434# Refactored over here to allow access by GSASIIscriptable.py
435def WriteAtomsNuclear(fp, phasedict, phasenam, parmDict, sigDict, labellist,
436                          RBparms={}):
437    'Write atom positions to CIF'
438    # phasedict = self.Phases[phasenam] # pointer to current phase info
439    General = phasedict['General']
440    cx,ct,cs,cia = General['AtomPtrs']
441    GS = G2lat.cell2GS(General['Cell'][1:7])
442    Amat = G2lat.cell2AB(General['Cell'][1:7])[0]
443    Atoms = phasedict['Atoms']
444    cfrac = cx+3
445    fpfx = str(phasedict['pId'])+'::Afrac:'
446    for i,at in enumerate(Atoms):
447        fval = parmDict.get(fpfx+str(i),at[cfrac])
448        if fval != 0.0:
449            break
450    else:
451        WriteCIFitem(fp, '\n# PHASE HAS NO ATOMS!')
452        return
453
454    WriteCIFitem(fp, '\n# ATOMIC COORDINATES AND DISPLACEMENT PARAMETERS')
455    WriteCIFitem(fp, 'loop_ '+
456                 '\n   _atom_site_label'+
457                 '\n   _atom_site_type_symbol'+
458                 '\n   _atom_site_fract_x'+
459                 '\n   _atom_site_fract_y'+
460                 '\n   _atom_site_fract_z'+
461                 '\n   _atom_site_occupancy'+
462                 '\n   _atom_site_adp_type'+
463                 '\n   _atom_site_U_iso_or_equiv'+
464                 '\n   _atom_site_site_symmetry_multiplicity')
465
466    varnames = {cx:'Ax',cx+1:'Ay',cx+2:'Az',cx+3:'Afrac',
467                cia+1:'AUiso',cia+2:'AU11',cia+3:'AU22',cia+4:'AU33',
468                cia+5:'AU12',cia+6:'AU13',cia+7:'AU23'}
469    # Empty the labellist
470    while labellist:
471        labellist.pop()
472
473    pfx = str(phasedict['pId'])+'::'
474    # loop over all atoms
475    naniso = 0
476    for i,at in enumerate(Atoms):
477        if phasedict['General']['Type'] == 'macromolecular':
478            label = '%s_%s_%s_%s'%(at[ct-1],at[ct-3],at[ct-4],at[ct-2])
479            s = PutInCol(MakeUniqueLabel(label,labellist),15) # label
480        else:
481            s = PutInCol(MakeUniqueLabel(at[ct-1],labellist),6) # label
482        fval = parmDict.get(fpfx+str(i),at[cfrac])
483        if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
484        s += PutInCol(FmtAtomType(at[ct]),4) # type
485        if at[cia] == 'I':
486            adp = 'Uiso '
487        else:
488            adp = 'Uani '
489            naniso += 1
490            t = G2lat.Uij2Ueqv(at[cia+2:cia+8],GS,Amat)[0]
491            for j in (2,3,4):
492                var = pfx+varnames[cia+j]+":"+str(i)
493        for j in (cx,cx+1,cx+2,cx+3,cia,cia+1):
494            if j in (cx,cx+1,cx+2):
495                dig = 11
496                sigdig = -0.00009
497            else:
498                dig = 10
499                sigdig = -0.0009
500            if j == cia:
501                s += adp
502            else:
503                var = pfx+varnames[j]+":"+str(i)
504                dvar = pfx+"d"+varnames[j]+":"+str(i)
505                if dvar not in sigDict:
506                    dvar = var
507                if j == cia+1 and adp == 'Uani ':
508                    sig = sigdig
509                    val = t
510                else:
511                    #print var,(var in parmDict),(var in sigDict)
512                    val = parmDict.get(var,at[j])
513                    sig = sigDict.get(dvar,sigdig)
514                    if dvar in G2mv.GetDependentVars(): # do not include an esd for dependent vars
515                        sig = -abs(sig)
516                s += PutInCol(G2mth.ValEsd(val,sig),dig)
517        s += PutInCol(at[cs+1],3)
518        WriteCIFitem(fp, s)
519    if naniso != 0: 
520        # now loop over aniso atoms
521        WriteCIFitem(fp, '\nloop_' + '\n   _atom_site_aniso_label' +
522                     '\n   _atom_site_aniso_U_11' + '\n   _atom_site_aniso_U_22' +
523                     '\n   _atom_site_aniso_U_33' + '\n   _atom_site_aniso_U_12' +
524                     '\n   _atom_site_aniso_U_13' + '\n   _atom_site_aniso_U_23')
525        for i,at in enumerate(Atoms):
526            fval = parmDict.get(fpfx+str(i),at[cfrac])
527            if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
528            if at[cia] == 'I': continue
529            s = PutInCol(labellist[i],6) # label
530            for j in (2,3,4,5,6,7):
531                sigdig = -0.0009
532                var = pfx+varnames[cia+j]+":"+str(i)
533                val = parmDict.get(var,at[cia+j])
534                sig = sigDict.get(var,sigdig)
535                s += PutInCol(G2mth.ValEsd(val,sig),11)
536            WriteCIFitem(fp, s)
537    # save information about rigid bodies
538    header = False
539    num = 0
540    rbAtoms = []
541    for irb,RBObj in enumerate(phasedict['RBModels'].get('Residue',[])):
542        if not header:
543            header = True
544            RBheader(fp)
545        jrb = RBparms['RBIds']['Residue'].index(RBObj['RBId'])
546        rbsx = str(irb)+':'+str(jrb)
547        num += 1
548        WriteCIFitem(fp,'',str(num))
549        RBModel = RBparms['Residue'][RBObj['RBId']]
550        SGData = phasedict['General']['SGData']
551        Sytsym,Mult = G2spc.SytSym(RBObj['Orig'][0],SGData)[:2]
552        s = '''GSAS-II residue rigid body "{}" with {} atoms
553  Site symmetry @ origin: {}, multiplicity: {}
554'''.format(RBObj['RBname'],len(RBModel['rbTypes']),Sytsym,Mult)
555        for i in G2stIO.WriteResRBModel(RBModel):
556            s += i
557        s += '\n Location:\n'
558        for i in G2stIO.WriteRBObjPOAndSig(pfx,'RBR',rbsx,parmDict,sigDict):
559            s += i+'\n'
560        for i in G2stIO.WriteRBObjTLSAndSig(pfx,'RBR',rbsx,
561                        RBObj['ThermalMotion'][0],parmDict,sigDict):
562            s += i
563        nTors = len(RBObj['Torsions'])
564        if nTors:
565            for i in G2stIO.WriteRBObjTorAndSig(pfx,rbsx,parmDict,sigDict,
566                        nTors):
567                s += i
568        WriteCIFitem(fp,'',s.rstrip())
569       
570        pId = phasedict['pId']
571        for i in RBObj['Ids']:
572            lbl = G2obj.LookupAtomLabel(pId,G2obj.LookupAtomId(pId,i))[0]
573            rbAtoms.append('{:7s} 1_555 {:3d} ?'.format(lbl,num))
574        #GSASIIpath.IPyBreak()
575
576    for irb,RBObj in enumerate(phasedict['RBModels'].get('Vector',[])):
577        if not header:
578            header = True
579            RBheader(fp)
580        jrb = RBparms['RBIds']['Vector'].index(RBObj['RBId'])
581        rbsx = str(irb)+':'+str(jrb)
582        num += 1
583        WriteCIFitem(fp,'',str(num))
584        RBModel = RBparms['Vector'][RBObj['RBId']]
585        SGData = phasedict['General']['SGData']
586        Sytsym,Mult = G2spc.SytSym(RBObj['Orig'][0],SGData)[:2]
587        s = '''GSAS-II vector rigid body "{}" with {} atoms
588  Site symmetry @ origin: {}, multiplicity: {}
589'''.format(RBObj['RBname'],len(RBModel['rbTypes']),Sytsym,Mult)
590        for i in G2stIO.WriteVecRBModel(RBModel,sigDict,irb):
591            s += i
592        s += '\n Location:\n'
593        for i in G2stIO.WriteRBObjPOAndSig(pfx,'RBV',rbsx,parmDict,sigDict):
594            s += i+'\n'
595        for i in G2stIO.WriteRBObjTLSAndSig(pfx,'RBV',rbsx,
596                        RBObj['ThermalMotion'][0],parmDict,sigDict):
597            s += i
598        WriteCIFitem(fp,'',s.rstrip())
599       
600        pId = phasedict['pId']
601        for i in RBObj['Ids']:
602            lbl = G2obj.LookupAtomLabel(pId,G2obj.LookupAtomId(pId,i))[0]
603            rbAtoms.append('{:7s} 1_555 {:3d} ?'.format(lbl,num))
604
605    if rbAtoms:
606        WriteCIFitem(fp,'loop_\n    _restr_rigid_body.id'+
607            '\n    _restr_rigid_body.atom_site_label\n    _restr_rigid_body.site_symmetry'+
608            '\n    _restr_rigid_body.class_id\n    _restr_rigid_body.details')
609        for i,l in enumerate(rbAtoms):
610            WriteCIFitem(fp,'   {:5d} {}'.format(i+1,l))
611           
612def WriteAtomsMagnetic(fp, phasedict, phasenam, parmDict, sigDict, labellist):
613    'Write atom positions to CIF'
614    # phasedict = self.Phases[phasenam] # pointer to current phase info
615    General = phasedict['General']
616    cx,ct,cs,cia = General['AtomPtrs']
617    Atoms = phasedict['Atoms']
618    cfrac = cx+3
619    fpfx = str(phasedict['pId'])+'::Afrac:'
620    for i,at in enumerate(Atoms):
621        fval = parmDict.get(fpfx+str(i),at[cfrac])
622        if fval != 0.0:
623            break
624    else:
625        WriteCIFitem(fp, '\n# PHASE HAS NO ATOMS!')
626        return
627
628    WriteCIFitem(fp, '\n# ATOMIC COORDINATES AND DISPLACEMENT PARAMETERS')
629    WriteCIFitem(fp, 'loop_ '+
630                 '\n   _atom_site_label'+
631                 '\n   _atom_site_type_symbol'+
632                 '\n   _atom_site_fract_x'+
633                 '\n   _atom_site_fract_y'+
634                 '\n   _atom_site_fract_z'+
635                 '\n   _atom_site_occupancy'+
636                 '\n   _atom_site_adp_type'+
637                 '\n   _atom_site_U_iso_or_equiv'+
638                 '\n   _atom_site_site_symmetry_multiplicity')
639
640    varnames = {cx:'Ax',cx+1:'Ay',cx+2:'Az',cx+3:'Afrac',
641                cx+4:'AMx',cx+5:'AMy',cx+6:'AMz',
642                cia+1:'AUiso',cia+2:'AU11',cia+3:'AU22',cia+4:'AU33',
643                cia+5:'AU12',cia+6:'AU13',cia+7:'AU23'}
644    # Empty the labellist
645    while labellist:
646        labellist.pop()
647
648    pfx = str(phasedict['pId'])+'::'
649    # loop over all atoms
650    naniso = 0
651    for i,at in enumerate(Atoms):
652        if phasedict['General']['Type'] == 'macromolecular':
653            label = '%s_%s_%s_%s'%(at[ct-1],at[ct-3],at[ct-4],at[ct-2])
654            s = PutInCol(MakeUniqueLabel(label,labellist),15) # label
655        else:
656            s = PutInCol(MakeUniqueLabel(at[ct-1],labellist),6) # label
657        fval = parmDict.get(fpfx+str(i),at[cfrac])
658        if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
659        s += PutInCol(FmtAtomType(at[ct]),4) # type
660        if at[cia] == 'I':
661            adp = 'Uiso '
662        else:
663            adp = 'Uani '
664            naniso += 1
665            # compute Uequiv crudely
666            # correct: Defined as "1/3 trace of diagonalized U matrix".
667            # SEE cell2GS & Uij2Ueqv to GSASIIlattice. Former is needed to make the GS matrix used by the latter.
668            t = 0.0
669            for j in (2,3,4):
670                var = pfx+varnames[cia+j]+":"+str(i)
671                t += parmDict.get(var,at[cia+j])
672        for j in (cx,cx+1,cx+2,cx+3,cia,cia+1):
673            if j in (cx,cx+1,cx+2):
674                dig = 11
675                sigdig = -0.00009
676            else:
677                dig = 10
678                sigdig = -0.009
679            if j == cia:
680                s += adp
681            else:
682                var = pfx+varnames[j]+":"+str(i)
683                dvar = pfx+"d"+varnames[j]+":"+str(i)
684                if dvar not in sigDict:
685                    dvar = var
686                if j == cia+1 and adp == 'Uani ':
687                    val = t/3.
688                    sig = sigdig
689                else:
690                    #print var,(var in parmDict),(var in sigDict)
691                    val = parmDict.get(var,at[j])
692                    sig = sigDict.get(dvar,sigdig)
693                    if dvar in G2mv.GetDependentVars(): # do not include an esd for dependent vars
694                        sig = -abs(sig)
695                s += PutInCol(G2mth.ValEsd(val,sig),dig)
696        s += PutInCol(at[cs+1],3)
697        WriteCIFitem(fp, s)
698    if naniso: 
699        # now loop over aniso atoms
700        WriteCIFitem(fp, '\nloop_' + '\n   _atom_site_aniso_label' +
701                     '\n   _atom_site_aniso_U_11' + '\n   _atom_site_aniso_U_22' +
702                     '\n   _atom_site_aniso_U_33' + '\n   _atom_site_aniso_U_12' +
703                     '\n   _atom_site_aniso_U_13' + '\n   _atom_site_aniso_U_23')
704        for i,at in enumerate(Atoms):
705            fval = parmDict.get(fpfx+str(i),at[cfrac])
706            if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
707            if at[cia] == 'I': continue
708            s = PutInCol(labellist[i],6) # label
709            for j in (2,3,4,5,6,7):
710                sigdig = -0.0009
711                var = pfx+varnames[cia+j]+":"+str(i)
712                val = parmDict.get(var,at[cia+j])
713                sig = sigDict.get(var,sigdig)
714                s += PutInCol(G2mth.ValEsd(val,sig),11)
715            WriteCIFitem(fp, s)
716    # now loop over mag atoms (e.g. all of them)
717    WriteCIFitem(fp, '\nloop_' + '\n   _atom_site_moment.label' +
718                 '\n   _atom_site_moment.crystalaxis_x' +
719                 '\n   _atom_site_moment.crystalaxis_y' +
720                 '\n   _atom_site_moment.crystalaxis_z')
721    for i,at in enumerate(Atoms):
722        fval = parmDict.get(fpfx+str(i),at[cfrac])
723        if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
724        s = PutInCol(labellist[i],6) # label
725        for j in (cx+4,cx+5,cx+6):
726            sigdig = -0.0009
727            var = pfx+varnames[j]+":"+str(i)
728            val = parmDict.get(var,at[j])
729            sig = sigDict.get(var,sigdig)
730            s += PutInCol(G2mth.ValEsd(val,sig),11)
731        WriteCIFitem(fp, s)
732
733def WriteAtomsMM(fp, phasedict, phasenam, parmDict, sigDict,
734                          RBparms={}):
735    'Write atom positions to CIF using mmCIF items'
736    AA3letter = ['ALA','ARG','ASN','ASP','CYS','GLN','GLU','GLY','HIS','ILE',
737                'LEU','LYS','MET','PHE','PRO','SER','THR','TRP','TYR','VAL','MSE']
738    # phasedict = self.Phases[phasenam] # pointer to current phase info
739    General = phasedict['General']
740    cx,ct,cs,cia = General['AtomPtrs']
741    #GS = G2lat.cell2GS(General['Cell'][1:7])
742    Amat = G2lat.cell2AB(General['Cell'][1:7])[0]
743    Atoms = phasedict['Atoms']
744    #cfrac = cx+3
745    #fpfx = str(phasedict['pId'])+'::Afrac:'
746    if len(Atoms) == 0:
747        WriteCIFitem(fp, '\n# PHASE HAS NO ATOMS!')
748        return
749
750    WriteCIFitem(fp, '\n# ATOMIC COORDINATES AND DISPLACEMENT PARAMETERS')
751    WriteCIFitem(fp, 'loop_ '+
752                 '\n   _atom_site.group_PDB'+
753                 '\n   _atom_site.id'+
754                 '\n   _atom_site.type_symbol'+
755                 '\n   _atom_site.label_atom_id'+
756                 '\n   _atom_site.auth_atom_id'+
757                 '\n   _atom_site.label_alt_id'+
758                 '\n   _atom_site.label_comp_id'+
759                 '\n   _atom_site.auth_comp_id'+
760                 '\n   _atom_site.label_asym_id'+
761                 '\n   _atom_site.auth_asym_id'+
762                 '\n   _atom_site.label_entity_id'+
763                 '\n   _atom_site.label_seq_id'+
764                 '\n   _atom_site.auth_seq_id'+
765                 '\n   _atom_site.pdbx_PDB_ins_code'+
766                 '\n   _atom_site.pdbx_formal_charge'+
767                 '\n   _atom_site.pdbx_PDB_model_num'
768                 '\n   _atom_site.fract_x'+
769                 '\n   _atom_site.fract_y'+
770                 '\n   _atom_site.fract_z'+
771                 '\n   _atom_site.occupancy'+
772                 '\n   _atom_site.B_iso_or_equiv'+
773                 '\n   _atom_site.Cartn_x'+
774                 '\n   _atom_site.Cartn_y'+
775                 '\n   _atom_site.Cartn_z'
776                 )
777
778# _atom_site.Cartn_x_esd
779# _atom_site.Cartn_y_esd
780# _atom_site.Cartn_z_esd
781# _atom_site.occupancy_esd
782#
783
784    varnames = {cx:'Ax',cx+1:'Ay',cx+2:'Az',cx+3:'Afrac',
785                cia+1:'AUiso',cia+2:'AU11',cia+3:'AU22',cia+4:'AU33',
786                cia+5:'AU12',cia+6:'AU13',cia+7:'AU23'}
787
788    pfx = str(phasedict['pId'])+'::'
789    num = 0
790    # uniquely index the side chains
791    entity_id = {}
792    for i,at in enumerate(Atoms): 
793        if at[ct-2] not in entity_id:
794            num += 1
795            entity_id[at[ct-2]] = num
796
797    # loop over all atoms
798#    naniso = 0
799    for i,at in enumerate(Atoms):
800        if at[ct-3] in AA3letter:
801            s = 'ATOM   '
802        else:
803            s = 'HETATM '
804        s += PutInCol(str(i+1),5) # atom number
805        s += PutInCol(FmtAtomType(at[ct]),4) # type
806        s += PutInCol(at[ct-1],4) # _atom_id
807        s += PutInCol(at[ct-1],4) # _atom_id
808        s += PutInCol('.',2) # alt_id
809        s += PutInCol(at[ct-3],4) # comp_id
810        s += PutInCol(at[ct-3],4) # comp_id
811        s += PutInCol(at[ct-2],3) # _asym_id
812        s += PutInCol(at[ct-2],3) # _asym_id
813        s += PutInCol(str(entity_id[at[ct-2]]),3) # entity_id
814        s += PutInCol(at[ct-4],2) # _seq_id
815        s += PutInCol(at[ct-4],2) # _seq_id
816        s += PutInCol('?',2) # pdbx_PDB_ins_code
817        s += PutInCol('?',2) # pdbx_formal_charge
818        s += PutInCol('1',2) # pdbx_PDB_model_num
819           
820#        fval = parmDict.get(fpfx+str(i),at[cfrac])
821#        if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
822#        if at[cia] == 'I':
823#            adp = 'Uiso '
824#        else:
825#            adp = 'Uani '
826#            naniso += 1
827#            t = G2lat.Uij2Ueqv(at[cia+2:cia+8],GS,Amat)[0]
828#            for j in (2,3,4):
829#                var = pfx+varnames[cia+j]+":"+str(i)
830        for j in (cx,cx+1,cx+2,cx+3,cia+1):
831            if j in (cx,cx+1,cx+2):
832                dig = 11
833                sigdig = -0.00009
834            else:
835                dig = 5
836                sigdig = -0.009
837            var = pfx+varnames[j]+":"+str(i)
838            dvar = pfx+"d"+varnames[j]+":"+str(i)
839            if dvar not in sigDict:
840                dvar = var
841            #print var,(var in parmDict),(var in sigDict)
842            val = parmDict.get(var,at[j])
843            sig = sigDict.get(dvar,sigdig)
844            if j == cia+1:  # convert U to B
845                val *= 8*np.pi**2
846                sig *= 8*np.pi**2
847            if dvar in G2mv.GetDependentVars(): # do not include an esd for dependent vars
848                sig = -abs(sig)
849            s += PutInCol(G2mth.ValEsd(val,sig),dig)
850        # Cartesian coordinates
851        for xyz in np.inner(Amat,at[cx:cx+3]):
852            s += PutInCol(G2mth.ValEsd(xyz,-0.009),8)
853        WriteCIFitem(fp, s)
854    # save information about rigid bodies
855#     header = False
856#     num = 0
857#     rbAtoms = []
858#     for irb,RBObj in enumerate(phasedict['RBModels'].get('Residue',[])):
859#         if not header:
860#             header = True
861#             RBheader(fp)
862#         jrb = RBparms['RBIds']['Residue'].index(RBObj['RBId'])
863#         rbsx = str(irb)+':'+str(jrb)
864#         num += 1
865#         WriteCIFitem(fp,'',str(num))
866#         RBModel = RBparms['Residue'][RBObj['RBId']]
867#         SGData = phasedict['General']['SGData']
868#         Sytsym,Mult = G2spc.SytSym(RBObj['Orig'][0],SGData)[:2]
869#         s = '''GSAS-II residue rigid body "{}" with {} atoms
870#   Site symmetry @ origin: {}, multiplicity: {}
871# '''.format(RBObj['RBname'],len(RBModel['rbTypes']),Sytsym,Mult)
872#         for i in G2stIO.WriteResRBModel(RBModel):
873#             s += i
874#         s += '\n Location:\n'
875#         for i in G2stIO.WriteRBObjPOAndSig(pfx,'RBR',rbsx,parmDict,sigDict):
876#             s += i+'\n'
877#         for i in G2stIO.WriteRBObjTLSAndSig(pfx,'RBR',rbsx,
878#                         RBObj['ThermalMotion'][0],parmDict,sigDict):
879#             s += i
880#         nTors = len(RBObj['Torsions'])
881#         if nTors:
882#             for i in G2stIO.WriteRBObjTorAndSig(pfx,rbsx,parmDict,sigDict,
883#                         nTors):
884#                 s += i
885#         WriteCIFitem(fp,'',s.rstrip())
886       
887#         pId = phasedict['pId']
888#         for i in RBObj['Ids']:
889#             lbl = G2obj.LookupAtomLabel(pId,G2obj.LookupAtomId(pId,i))[0]
890#             rbAtoms.append('{:7s} 1_555 {:3d} ?'.format(lbl,num))
891#         #GSASIIpath.IPyBreak()
892
893#     for irb,RBObj in enumerate(phasedict['RBModels'].get('Vector',[])):
894#         if not header:
895#             header = True
896#             RBheader(fp)
897#         jrb = RBparms['RBIds']['Vector'].index(RBObj['RBId'])
898#         rbsx = str(irb)+':'+str(jrb)
899#         num += 1
900#         WriteCIFitem(fp,'',str(num))
901#         RBModel = RBparms['Vector'][RBObj['RBId']]
902#         SGData = phasedict['General']['SGData']
903#         Sytsym,Mult = G2spc.SytSym(RBObj['Orig'][0],SGData)[:2]
904#         s = '''GSAS-II vector rigid body "{}" with {} atoms
905#   Site symmetry @ origin: {}, multiplicity: {}
906# '''.format(RBObj['RBname'],len(RBModel['rbTypes']),Sytsym,Mult)
907#         for i in G2stIO.WriteVecRBModel(RBModel,sigDict,irb):
908#             s += i
909#         s += '\n Location:\n'
910#         for i in G2stIO.WriteRBObjPOAndSig(pfx,'RBV',rbsx,parmDict,sigDict):
911#             s += i+'\n'
912#         for i in G2stIO.WriteRBObjTLSAndSig(pfx,'RBV',rbsx,
913#                         RBObj['ThermalMotion'][0],parmDict,sigDict):
914#             s += i
915#         WriteCIFitem(fp,'',s.rstrip())
916       
917#         pId = phasedict['pId']
918#         for i in RBObj['Ids']:
919#             lbl = G2obj.LookupAtomLabel(pId,G2obj.LookupAtomId(pId,i))[0]
920#             rbAtoms.append('{:7s} 1_555 {:3d} ?'.format(lbl,num))
921
922#     if rbAtoms:
923#         WriteCIFitem(fp,'loop_\n    _restr_rigid_body.id'+
924#             '\n    _restr_rigid_body.atom_site_label\n    _restr_rigid_body.site_symmetry'+
925#             '\n    _restr_rigid_body.class_id\n    _restr_rigid_body.details')
926#         for i,l in enumerate(rbAtoms):
927#             WriteCIFitem(fp,'   {:5d} {}'.format(i+1,l))
928
929# Refactored over here to allow access by GSASIIscriptable.py
930def WriteSeqAtomsNuclear(fp, cell, phasedict, phasenam, hist, seqData, RBparms):
931    'Write atom positions to CIF'
932    General = phasedict['General']
933    cx,ct,cs,cia = General['AtomPtrs']
934    GS = G2lat.cell2GS(cell[:6])
935    Amat = G2lat.cell2AB(cell[:6])[0]
936
937    # phasedict = self.Phases[phasenam] # pointer to current phase info
938    parmDict = seqData[hist]['parmDict']
939    sigDict = dict(zip(seqData[hist]['varyList'],seqData[hist]['sig']))
940    Atoms = phasedict['Atoms']
941    cfrac = cx+3
942    fpfx = str(phasedict['pId'])+'::Afrac:'
943    for i,at in enumerate(Atoms):
944        fval = parmDict.get(fpfx+str(i),at[cfrac])
945        if fval != 0.0:
946            break
947    else:
948        WriteCIFitem(fp, '\n# PHASE HAS NO ATOMS!')
949        return
950
951    WriteCIFitem(fp, '\n# ATOMIC COORDINATES AND DISPLACEMENT PARAMETERS')
952    WriteCIFitem(fp, 'loop_ '+
953                 '\n   _atom_site_label'+
954                 '\n   _atom_site_type_symbol'+
955                 '\n   _atom_site_fract_x'+
956                 '\n   _atom_site_fract_y'+
957                 '\n   _atom_site_fract_z'+
958                 '\n   _atom_site_occupancy'+
959                 '\n   _atom_site_adp_type'+
960                 '\n   _atom_site_U_iso_or_equiv'+
961                 '\n   _atom_site_site_symmetry_multiplicity')
962
963    varnames = {cx:'Ax',cx+1:'Ay',cx+2:'Az',cx+3:'Afrac',
964                cia+1:'AUiso',cia+2:'AU11',cia+3:'AU22',cia+4:'AU33',
965                cia+5:'AU12',cia+6:'AU13',cia+7:'AU23'}
966
967    labellist = []  # used to make atom labels unique as required in CIF
968    pfx = str(phasedict['pId'])+'::'
969    # loop over all atoms
970    naniso = 0
971
972    for i,at in enumerate(Atoms):
973        if phasedict['General']['Type'] == 'macromolecular':
974            label = '%s_%s_%s_%s'%(at[ct-1],at[ct-3],at[ct-4],at[ct-2])
975            s = PutInCol(MakeUniqueLabel(label,labellist),15) # label
976        else:
977            s = PutInCol(MakeUniqueLabel(at[ct-1],labellist),6) # label
978        fval = parmDict.get(fpfx+str(i),at[cfrac])
979        if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
980        s += PutInCol(FmtAtomType(at[ct]),4) # type
981        if at[cia] == 'I':
982            adp = 'Uiso '
983        else:
984            adp = 'Uani '
985            naniso += 1
986            t = G2lat.Uij2Ueqv(at[cia+2:cia+8],GS,Amat)[0]
987            for j in (2,3,4):
988                var = pfx+varnames[cia+j]+":"+str(i)
989        for j in (cx,cx+1,cx+2,cx+3,cia,cia+1):
990            if j in (cx,cx+1,cx+2):
991                dig = 11
992                sigdig = -0.00009
993            else:
994                dig = 10
995                sigdig = -0.0009
996            if j == cia:
997                s += adp
998            else:
999                var = pfx+varnames[j]+":"+str(i)
1000                dvar = pfx+"d"+varnames[j]+":"+str(i)
1001                if dvar not in sigDict:
1002                    dvar = var
1003                if j == cia+1 and adp == 'Uani ':
1004                    sig = sigdig
1005                    val = t
1006                else:
1007                    #print var,(var in parmDict),(var in sigDict)
1008                    val = parmDict.get(var,at[j])
1009                    sig = sigDict.get(dvar,sigdig)
1010                    if dvar in G2mv.GetDependentVars(): # do not include an esd for dependent vars
1011                        sig = -abs(sig)
1012                s += PutInCol(G2mth.ValEsd(val,sig),dig)
1013        s += PutInCol(at[cs+1],3)
1014        WriteCIFitem(fp, s)
1015    if naniso != 0: 
1016        # now loop over aniso atoms
1017        WriteCIFitem(fp, '\nloop_' + '\n   _atom_site_aniso_label' +
1018                     '\n   _atom_site_aniso_U_11' + '\n   _atom_site_aniso_U_22' +
1019                     '\n   _atom_site_aniso_U_33' + '\n   _atom_site_aniso_U_12' +
1020                     '\n   _atom_site_aniso_U_13' + '\n   _atom_site_aniso_U_23')
1021        for i,at in enumerate(Atoms):
1022            fval = parmDict.get(fpfx+str(i),at[cfrac])
1023            if fval == 0.0: continue # ignore any atoms that have a occupancy set to 0 (exact)
1024            if at[cia] == 'I': continue
1025            s = PutInCol(labellist[i],6) # label
1026            for j in (2,3,4,5,6,7):
1027                sigdig = -0.0009
1028                var = pfx+varnames[cia+j]+":"+str(i)
1029                val = parmDict.get(var,at[cia+j])
1030                sig = sigDict.get(var,sigdig)
1031                s += PutInCol(G2mth.ValEsd(val,sig),11)
1032            WriteCIFitem(fp, s)
1033    # save information about rigid bodies
1034    header = False
1035    num = 0
1036    rbAtoms = []
1037    for irb,RBObj in enumerate(phasedict['RBModels'].get('Residue',[])):
1038        if not header:
1039            header = True
1040            RBheader(fp)
1041        jrb = RBparms['RBIds']['Residue'].index(RBObj['RBId'])
1042        rbsx = str(irb)+':'+str(jrb)
1043        num += 1
1044        WriteCIFitem(fp,'',str(num))
1045        RBModel = RBparms['Residue'][RBObj['RBId']]
1046        SGData = phasedict['General']['SGData']
1047        Sytsym,Mult = G2spc.SytSym(RBObj['Orig'][0],SGData)[:2]
1048        s = '''GSAS-II residue rigid body "{}" with {} atoms
1049  Site symmetry @ origin: {}, multiplicity: {}
1050'''.format(RBObj['RBname'],len(RBModel['rbTypes']),Sytsym,Mult)
1051        for i in G2stIO.WriteResRBModel(RBModel):
1052            s += i
1053        s += '\n Location:\n'
1054        for i in G2stIO.WriteRBObjPOAndSig(pfx,'RBR',rbsx,parmDict,sigDict):
1055            s += i+'\n'
1056        for i in G2stIO.WriteRBObjTLSAndSig(pfx,'RBR',rbsx,
1057                        RBObj['ThermalMotion'][0],parmDict,sigDict):
1058            s += i
1059        nTors = len(RBObj['Torsions'])
1060        if nTors:
1061            for i in G2stIO.WriteRBObjTorAndSig(pfx,rbsx,parmDict,sigDict,
1062                        nTors):
1063                s += i
1064        WriteCIFitem(fp,'',s.rstrip())
1065       
1066        pId = phasedict['pId']
1067        for i in RBObj['Ids']:
1068            lbl = G2obj.LookupAtomLabel(pId,G2obj.LookupAtomId(pId,i))[0]
1069            rbAtoms.append('{:7s} 1_555 {:3d} ?'.format(lbl,num))
1070        #GSASIIpath.IPyBreak()
1071
1072    for irb,RBObj in enumerate(phasedict['RBModels'].get('Vector',[])):
1073        if not header:
1074            header = True
1075            RBheader(fp)
1076        jrb = RBparms['RBIds']['Vector'].index(RBObj['RBId'])
1077        rbsx = str(irb)+':'+str(jrb)
1078        num += 1
1079        WriteCIFitem(fp,'',str(num))
1080        RBModel = RBparms['Vector'][RBObj['RBId']]
1081        SGData = phasedict['General']['SGData']
1082        Sytsym,Mult = G2spc.SytSym(RBObj['Orig'][0],SGData)[:2]
1083        s = '''GSAS-II vector rigid body "{}" with {} atoms
1084  Site symmetry @ origin: {}, multiplicity: {}
1085'''.format(RBObj['RBname'],len(RBModel['rbTypes']),Sytsym,Mult)
1086        for i in G2stIO.WriteVecRBModel(RBModel,sigDict,irb):
1087            s += i
1088        s += '\n Location:\n'
1089        for i in G2stIO.WriteRBObjPOAndSig(pfx,'RBV',rbsx,parmDict,sigDict):
1090            s += i+'\n'
1091        for i in G2stIO.WriteRBObjTLSAndSig(pfx,'RBV',rbsx,
1092                        RBObj['ThermalMotion'][0],parmDict,sigDict):
1093            s += i
1094        WriteCIFitem(fp,'',s.rstrip())
1095       
1096        pId = phasedict['pId']
1097        for i in RBObj['Ids']:
1098            lbl = G2obj.LookupAtomLabel(pId,G2obj.LookupAtomId(pId,i))[0]
1099            rbAtoms.append('{:7s} 1_555 {:3d} ?'.format(lbl,num))
1100
1101    if rbAtoms:
1102        WriteCIFitem(fp,'loop_\n    _restr_rigid_body.id'+
1103            '\n    _restr_rigid_body.atom_site_label\n    _restr_rigid_body.site_symmetry'+
1104            '\n    _restr_rigid_body.class_id\n    _restr_rigid_body.details')
1105        for i,l in enumerate(rbAtoms):
1106            WriteCIFitem(fp,'   {:5d} {}'.format(i+1,l))
1107           
1108# Refactored over here to allow access by GSASIIscriptable.py
1109def MakeUniqueLabel(lbl, labellist):
1110    lbl = lbl.strip()
1111    if not lbl: # deal with a blank label
1112        lbl = 'A_1'
1113    if lbl not in labellist:
1114        labellist.append(lbl)
1115        return lbl
1116    i = 1
1117    prefix = lbl
1118    if '_' in lbl:
1119        prefix = lbl[:lbl.rfind('_')]
1120        suffix = lbl[lbl.rfind('_')+1:]
1121        try:
1122            i = int(suffix)+1
1123        except:
1124            pass
1125    while prefix+'_'+str(i) in labellist:
1126        i += 1
1127    else:
1128        lbl = prefix+'_'+str(i)
1129        labellist.append(lbl)
1130
1131
1132# Refactored over here to allow access by GSASIIscriptable.py
1133def HillSortElements(elmlist):
1134    '''Sort elements in "Hill" order: C, H, others, (where others
1135    are alphabetical).
1136
1137    :params list elmlist: a list of element strings
1138
1139    :returns: a sorted list of element strings
1140    '''
1141    newlist = []
1142    oldlist = elmlist[:]
1143    for elm in ('C','H'):
1144        if elm in elmlist:
1145            newlist.append(elm)
1146            oldlist.pop(oldlist.index(elm))
1147    return newlist+sorted(oldlist)
1148
1149
1150def FmtAtomType(sym):
1151    'Reformat a GSAS-II atom type symbol to match CIF rules'
1152    sym = sym.replace('_','') # underscores are not allowed: no isotope designation?
1153    # in CIF, oxidation state sign symbols come after, not before
1154    if '+' in sym:
1155        sym = sym.replace('+','') + '+'
1156    elif '-' in sym:
1157        sym = sym.replace('-','') + '-'
1158    return sym
1159
1160
1161def PutInCol(val, wid):
1162    val = str(val).replace(' ', '')
1163    if not val: val = '?'
1164    fmt = '{:' + str(wid) + '} '
1165    try:
1166        return fmt.format(val)
1167    except TypeError:
1168        return fmt.format('.')
1169
1170
1171# Refactored over here to allow access by GSASIIscriptable.py
1172def WriteComposition(fp, phasedict, phasenam, parmDict, quickmode=True, keV=None):
1173    '''determine the composition for the unit cell, crudely determine Z and
1174    then compute the composition in formula units.
1175
1176    If quickmode is False, then scattering factors are added to the element loop.
1177
1178    If keV is specified, then resonant scattering factors are also computed and included.
1179    '''
1180    General = phasedict['General']
1181    Z = General.get('cellZ',0.0)
1182    cx,ct,cs,cia = General['AtomPtrs']
1183    Atoms = phasedict['Atoms']
1184    fpfx = str(phasedict['pId'])+'::Afrac:'
1185    cfrac = cx+3
1186    cmult = cs+1
1187    compDict = {} # combines H,D & T
1188    sitemultlist = []
1189    massDict = dict(zip(General['AtomTypes'],General['AtomMass']))
1190    cellmass = 0
1191    elmLookup = {}
1192    for i,at in enumerate(Atoms):
1193        atype = at[ct].strip()
1194        if atype.find('-') != -1: atype = atype.split('-')[0]
1195        if atype.find('+') != -1: atype = atype.split('+')[0]
1196        atype = atype[0].upper()+atype[1:2].lower() # force case conversion
1197        if atype == "D" or atype == "D": atype = "H"
1198        fvar = fpfx+str(i)
1199        fval = parmDict.get(fvar,at[cfrac])
1200        mult = at[cmult]
1201        if not massDict.get(at[ct]):
1202            print('Error: No mass found for atom type '+at[ct])
1203            print('Will not compute cell contents for phase '+phasenam)
1204            return
1205        cellmass += massDict[at[ct]]*mult*fval
1206        compDict[atype] = compDict.get(atype,0.0) + mult*fval
1207        elmLookup[atype] = at[ct].strip()
1208        if fval == 1: sitemultlist.append(mult)
1209    if len(compDict.keys()) == 0: return # no elements!
1210    if Z < 1: # Z has not been computed or set by user
1211        Z = 1
1212        if not sitemultlist:
1213            General['cellZ'] = 1
1214            return
1215        for i in range(2,min(sitemultlist)+1):
1216            for m in sitemultlist:
1217                if m % i != 0:
1218                    break
1219                else:
1220                    Z = i
1221        General['cellZ'] = Z # save it
1222
1223    if not quickmode:
1224        FFtable = G2el.GetFFtable(General['AtomTypes'])
1225        BLtable = G2el.GetBLtable(General)
1226
1227    WriteCIFitem(fp, '\nloop_  _atom_type_symbol _atom_type_number_in_cell')
1228    s = '       '
1229    if not quickmode:
1230        for j in ('a1','a2','a3','a4','b1','b2','b3','b4','c',2,1):
1231            if len(s) > 80:
1232                WriteCIFitem(fp, s)
1233                s = '       '
1234            if j==1:
1235                s += ' _atom_type_scat_source'
1236            elif j==2:
1237                s += ' _atom_type_scat_length_neutron'
1238            else:
1239                s += ' _atom_type_scat_Cromer_Mann_'
1240                s += j
1241        if keV:
1242            WriteCIFitem(fp, s)
1243            s = '        _atom_type_scat_dispersion_real _atom_type_scat_dispersion_imag _atom_type_scat_dispersion_source'
1244        WriteCIFitem(fp, s)
1245
1246           
1247    formula = ''
1248    for elem in HillSortElements(list(compDict.keys())):
1249        s = '  '
1250        elmsym = elmLookup[elem]
1251        # CIF does not allow underscore in element symbol (https://www.iucr.org/__data/iucr/cifdic_html/1/cif_core.dic/Iatom_type_symbol.html)
1252        if elmsym.endswith("_"):
1253             s += PutInCol(elmsym.replace('_',''))
1254        elif '_' in elmsym:
1255             s += PutInCol(elmsym.replace('_','~'))
1256        else:
1257            s += PutInCol(elmsym,7)
1258        s += PutInCol(G2mth.ValEsd(compDict[elem],-0.009,True),5)
1259        if not quickmode:
1260            for i in 'fa','fb','fc':
1261                if i != 'fc':
1262                    for j in range(4):
1263                        if elmsym in FFtable:
1264                            val = G2mth.ValEsd(FFtable[elmsym][i][j],-0.0009,True)
1265                        else:
1266                            val = '?'
1267                        s += ' '
1268                        s += PutInCol(val,9)
1269                else:
1270                    if elmsym in FFtable:
1271                        val = G2mth.ValEsd(FFtable[elmsym][i],-0.0009,True)
1272                    else:
1273                        val = '?'
1274                    s += ' '
1275                    s += PutInCol(val,9)
1276            if elmsym in BLtable:
1277                bldata = BLtable[elmsym]
1278                #isotope = bldata[0]
1279                #mass = bldata[1]['Mass']
1280                if 'BW-LS' in bldata[1]:
1281                    val = 0
1282                else:
1283                    val = G2mth.ValEsd(bldata[1]['SL'][0],-0.0009,True)
1284            else:
1285                val = '?'
1286            s += ' '
1287            s += PutInCol(val,9)
1288            WriteCIFitem(fp,s.rstrip())
1289            WriteCIFitem(fp,'      https://subversion.xray.aps.anl.gov/pyGSAS/trunk/atmdata.py')
1290            if keV:
1291                Orbs = G2el.GetXsectionCoeff(elem.split('+')[0].split('-')[0])
1292                FP,FPP,Mu = G2el.FPcalc(Orbs, keV)
1293                WriteCIFitem(fp,{:8.3f}{:8.3f}   https://subversion.xray.aps.anl.gov/pyGSAS/trunk/atmdata.py'.format(FP,FPP))
1294        else:
1295            WriteCIFitem(fp,s.rstrip())
1296        if formula: formula += " "
1297        formula += elem
1298        if compDict[elem] == Z: continue
1299        formula += G2mth.ValEsd(compDict[elem]/Z,-0.009,True)
1300    WriteCIFitem(fp,  '\n# Note that Z affects _cell_formula_sum and _weight')
1301    WriteCIFitem(fp,  '_cell_formula_units_Z',str(Z))
1302    WriteCIFitem(fp,  '_chemical_formula_sum',formula)
1303    WriteCIFitem(fp,  '_chemical_formula_weight',
1304                  G2mth.ValEsd(cellmass/Z,-0.09,True))
1305
1306def WriteCompositionMM(fp, phasedict, phasenam, parmDict, quickmode=True, keV=None):
1307    '''determine the composition for the unit cell, crudely determine Z and
1308    then compute the composition in formula units.
1309
1310    If quickmode is False, then scattering factors are added to the element loop.
1311
1312    If keV is specified, then resonant scattering factors are also computed and included.
1313    '''
1314    General = phasedict['General']
1315    Z = General.get('cellZ',0.0)
1316    cx,ct,cs,cia = General['AtomPtrs']
1317    Atoms = phasedict['Atoms']
1318    fpfx = str(phasedict['pId'])+'::Afrac:'
1319    cfrac = cx+3
1320    cmult = cs+1
1321    compDict = {} # combines H,D & T
1322    sitemultlist = []
1323    massDict = dict(zip(General['AtomTypes'],General['AtomMass']))
1324    cellmass = 0
1325    elmLookup = {}
1326    for i,at in enumerate(Atoms):
1327        atype = at[ct].strip()
1328        if atype.find('-') != -1: atype = atype.split('-')[0]
1329        if atype.find('+') != -1: atype = atype.split('+')[0]
1330        atype = atype[0].upper()+atype[1:2].lower() # force case conversion
1331        if atype == "D" or atype == "D": atype = "H"
1332        fvar = fpfx+str(i)
1333        fval = parmDict.get(fvar,at[cfrac])
1334        mult = at[cmult]
1335        if not massDict.get(at[ct]):
1336            print('Error: No mass found for atom type '+at[ct])
1337            print('Will not compute cell contents for phase '+phasenam)
1338            return
1339        cellmass += massDict[at[ct]]*mult*fval
1340        compDict[atype] = compDict.get(atype,0.0) + mult*fval
1341        elmLookup[atype] = at[ct].strip()
1342        if fval == 1: sitemultlist.append(mult)
1343    if len(compDict.keys()) == 0: return # no elements!
1344    if Z < 1: # Z has not been computed or set by user
1345        Z = 1
1346        if not sitemultlist:
1347            General['cellZ'] = 1
1348            return
1349        for i in range(2,min(sitemultlist)+1):
1350            for m in sitemultlist:
1351                if m % i != 0:
1352                    break
1353                else:
1354                    Z = i
1355        General['cellZ'] = Z # save it
1356
1357    if not quickmode:
1358        FFtable = G2el.GetFFtable(General['AtomTypes'])
1359        BLtable = G2el.GetBLtable(General)
1360
1361    WriteCIFitem(fp, '\nloop_  _atom_type.symbol _atom_type.number_in_cell')
1362    s = '       '
1363    if not quickmode:
1364        for j in ('a1','a2','a3','a4','b1','b2','b3','b4','c',2,1):
1365            if len(s) > 80:
1366                WriteCIFitem(fp, s)
1367                s = '       '
1368            if j==1:
1369                s += ' _atom_type.scat_source'
1370            elif j==2:
1371                s += ' _atom_type.scat_length_neutron'
1372            else:
1373                s += ' _atom_type.scat_Cromer_Mann_'
1374                s += j
1375        if keV:
1376            WriteCIFitem(fp, s)
1377            s = '        _atom_type.scat_dispersion_real _atom_type.scat_dispersion_imag _atom_type_scat_dispersion_source'
1378        WriteCIFitem(fp, s)
1379
1380           
1381    formula = ''
1382    for elem in HillSortElements(list(compDict.keys())):
1383        s = '  '
1384        elmsym = elmLookup[elem]
1385        # CIF does not allow underscore in element symbol (https://www.iucr.org/__data/iucr/cifdic_html/1/cif_core.dic/Iatom_type_symbol.html)
1386        if elmsym.endswith("_"):
1387             s += PutInCol(elmsym.replace('_',''))
1388        elif '_' in elmsym:
1389             s += PutInCol(elmsym.replace('_','~'))
1390        else:
1391            s += PutInCol(elmsym,7)
1392        s += PutInCol(G2mth.ValEsd(compDict[elem],-0.009,True),5)
1393        if not quickmode:
1394            for i in 'fa','fb','fc':
1395                if i != 'fc':
1396                    for j in range(4):
1397                        if elmsym in FFtable:
1398                            val = G2mth.ValEsd(FFtable[elmsym][i][j],-0.0009,True)
1399                        else:
1400                            val = '?'
1401                        s += ' '
1402                        s += PutInCol(val,9)
1403                else:
1404                    if elmsym in FFtable:
1405                        val = G2mth.ValEsd(FFtable[elmsym][i],-0.0009,True)
1406                    else:
1407                        val = '?'
1408                    s += ' '
1409                    s += PutInCol(val,9)
1410            if elmsym in BLtable:
1411                bldata = BLtable[elmsym]
1412                #isotope = bldata[0]
1413                #mass = bldata[1]['Mass']
1414                if 'BW-LS' in bldata[1]:
1415                    val = 0
1416                else:
1417                    val = G2mth.ValEsd(bldata[1]['SL'][0],-0.0009,True)
1418            else:
1419                val = '?'
1420            s += ' '
1421            s += PutInCol(val,9)
1422            WriteCIFitem(fp,s.rstrip())
1423            WriteCIFitem(fp,'      https://subversion.xray.aps.anl.gov/pyGSAS/trunk/atmdata.py')
1424            if keV:
1425                Orbs = G2el.GetXsectionCoeff(elem.split('+')[0].split('-')[0])
1426                FP,FPP,Mu = G2el.FPcalc(Orbs, keV)
1427                WriteCIFitem(fp,{:8.3f}{:8.3f}   https://subversion.xray.aps.anl.gov/pyGSAS/trunk/atmdata.py'.format(FP,FPP))
1428        else:
1429            WriteCIFitem(fp,s.rstrip())
1430        if formula: formula += " "
1431        formula += elem
1432        if compDict[elem] == Z: continue
1433        formula += G2mth.ValEsd(compDict[elem]/Z,-0.009,True)
1434    WriteCIFitem(fp,  '\n# Note that Z affects _cell_formula.sum and .weight')
1435    WriteCIFitem(fp,  '_cell.formula_units_Z',str(Z))
1436    WriteCIFitem(fp,  '_chemical_formula.sum',formula)
1437    WriteCIFitem(fp,  '_chemical_formula.weight',
1438                  G2mth.ValEsd(cellmass/Z,-0.09,True))
1439
1440class ExportCIF(G2IO.ExportBaseclass):
1441    '''Base class for CIF exports
1442    '''
1443    def __init__(self,G2frame,formatName,extension,longFormatName=None,):
1444        G2IO.ExportBaseclass.__init__(self,G2frame,formatName,extension,longFormatName=None)
1445        self.exporttype = []
1446        self.author = ''
1447        self.CIFname = ''
1448
1449    def ValidateAscii(self,checklist):
1450        '''Validate items as ASCII'''
1451        msg = ''
1452        for lbl,val in checklist:
1453            if not all(ord(c) < 128 for c in val):
1454                if msg: msg += '\n'
1455                msg += lbl + " contains unicode characters: " + val
1456        if msg:
1457            G2G.G2MessageBox(self.G2frame,
1458                             'Error: CIFs can contain only ASCII characters. Please change item(s) below:\n\n'+msg,
1459                             'Unicode not valid for CIF')
1460            return True
1461
1462    def _CellSelectNeeded(self,phasenam):
1463        '''Determines if selection is needed for a T value in a multiblock CIF
1464
1465        :returns: True if the choice of T is ambiguous and a human should
1466          be asked.
1467        '''
1468        phasedict = self.Phases[phasenam] # pointer to current phase info
1469        Tlist = {}  # histname & T values used for cell w/o Hstrain
1470        DijTlist = {} # hId & T values used for cell w/Hstrain
1471        # scan over histograms used in this phase to determine the best
1472        # data collection T value
1473        for h in phasedict['Histograms']:
1474            if not phasedict['Histograms'][h]['Use']: continue
1475            if 'Flack' in phasedict['Histograms'][h].keys():       #single crystal data
1476                return False
1477            T = self.Histograms[h]['Sample Parameters']['Temperature']
1478            if np.any(abs(np.array(phasedict['Histograms'][h]['HStrain'][0])) > 1e-8):
1479                DijTlist[h] = T
1480            else:
1481                Tlist[h] = T
1482        if len(Tlist) > 0:
1483            T = sum(Tlist.values())/len(Tlist)
1484            if max(Tlist.values()) - T > 1:
1485                return True # temperatures span more than 1 degree, user needs to pick one
1486            return False
1487        elif len(DijTlist) == 1:
1488            return False
1489        elif len(DijTlist) > 1:
1490                # each histogram has different cell lengths, user needs to pick one
1491            return True
1492       
1493    def _CellSelectHist(self,phasenam):
1494        '''Select T value for a phase in a multiblock CIF
1495
1496        :returns: T,h_ranId where T is a temperature (float) or '?' and
1497          h_ranId is the random Id (ranId) for a histogram in the
1498          current phase. This is stored in OverallParms['Controls']['CellHistSelection']
1499        '''
1500        phasedict = self.Phases[phasenam] # pointer to current phase info
1501        Tlist = {}  # histname & T values used for cell w/o Hstrain
1502        DijTlist = {} # hId & T values used for cell w/Hstrain
1503        # scan over histograms used in this phase to determine the best
1504        # data collection T value
1505        for h in phasedict['Histograms']:
1506            if not phasedict['Histograms'][h]['Use']: continue
1507            if 'Flack' in phasedict['Histograms'][h].keys():       #single crystal data
1508                return (300,None)
1509            T = self.Histograms[h]['Sample Parameters']['Temperature']
1510            if np.any(abs(np.array(phasedict['Histograms'][h]['HStrain'][0])) > 1e-8):
1511                DijTlist[h] = T
1512            else:
1513                Tlist[h] = T
1514        if len(Tlist) > 0:
1515            T = sum(Tlist.values())/len(Tlist)
1516            if max(Tlist.values()) - T > 1:
1517                # temperatures span more than 1 degree, user needs to pick one
1518                choices = ["{} (unweighted average)".format(T)]
1519                Ti = [T]
1520                for h in Tlist:
1521                    choices += ["{} (hist {})".format(Tlist[h],h)]
1522                    Ti += [Tlist[h]]                           
1523                msg = 'The cell parameters for phase {} are from\nhistograms with different temperatures.\n\nSelect a T value below'.format(phasenam)
1524                dlg = wx.SingleChoiceDialog(self.G2frame,msg,'Select T',choices)
1525                if dlg.ShowModal() == wx.ID_OK:
1526                    T = Ti[dlg.GetSelection()]
1527                else:
1528                    T = '?'
1529                dlg.Destroy()
1530            return (T,None)
1531        elif len(DijTlist) == 1:
1532            h = list(DijTlist.keys())[0]
1533            h_ranId = self.Histograms[h]['ranId']
1534            return (DijTlist[h],h_ranId)
1535        elif len(DijTlist) > 1:
1536            # each histogram has different cell lengths, user needs to pick one
1537            choices = []
1538            hi = []
1539            for h in DijTlist:
1540                choices += ["{} (hist {})".format(DijTlist[h],h)]
1541                hi += [h]
1542            msg = 'There are {} sets of cell parameters for phase {}\n due to refined Hstrain values.\n\nSelect the histogram to use with the phase form list below'.format(len(DijTlist),phasenam)
1543            dlg = wx.SingleChoiceDialog(self.G2frame,msg,'Select cell',choices)
1544            if dlg.ShowModal() == wx.ID_OK:
1545                h = hi[dlg.GetSelection()] 
1546                h_ranId = self.Histograms[h]['ranId']
1547                T = DijTlist[h]
1548            else:
1549                T = '?'
1550                h_ranId = None
1551            dlg.Destroy()
1552            return (T,h_ranId)
1553        else:
1554            print('Unexpected option in _CellSelectHist for',phasenam)
1555            return ('?',None)
1556
1557    def ShowHstrainCells(self,phasenam,datablockidDict):
1558        '''Displays the unit cell parameters for phases where Dij values create
1559        mutiple sets of lattice parameters. At present there is no way defined for this in
1560        CIF, so local data names are used.
1561        '''
1562        phasedict = self.Phases[phasenam] # pointer to current phase info
1563        Tlist = {}  # histname & T values used for cell w/o Hstrain
1564        DijTlist = {} # hId & T values used for cell w/Hstrain
1565        # scan over histograms used in this phase
1566        for h in phasedict['Histograms']:
1567            if not phasedict['Histograms'][h]['Use']: continue
1568            if np.any(abs(np.array(phasedict['Histograms'][h]['HStrain'][0])) > 1e-8):
1569                DijTlist[h] = self.Histograms[h]['Sample Parameters']['Temperature']
1570            else:
1571                Tlist[h] = self.Histograms[h]['Sample Parameters']['Temperature']
1572        if len(DijTlist) == 0: return
1573        if len(Tlist) + len(DijTlist) < 2: return
1574        SGData = phasedict['General']['SGData']
1575        for i in range(len(G2py3.cellGUIlist)):
1576            if SGData['SGLaue'] in G2py3.cellGUIlist[i][0]:
1577                terms = G2py3.cellGUIlist[i][5] + [6]
1578                break
1579        else:
1580            print('ShowHstrainCells error: Laue class not found',SGData['SGLaue'])
1581            terms = list(range(7))
1582       
1583        WriteCIFitem(self.fp, '\n# cell parameters generated by hydrostatic strain')
1584        WriteCIFitem(self.fp, 'loop_')
1585        WriteCIFitem(self.fp, '\t _gsas_measurement_temperature')
1586        for i in terms:
1587            WriteCIFitem(self.fp, '\t _gsas_cell_'+cellNames[i])
1588        WriteCIFitem(self.fp, '\t _gsas_cell_histogram_blockid')
1589        for h,T in Tlist.items():
1590            pId = phasedict['pId']
1591            hId = self.Histograms[h]['hId']
1592            cellList,cellSig = G2stIO.getCellSU(pId,hId,
1593                                        phasedict['General']['SGData'],
1594                                        self.parmDict,
1595                                        self.OverallParms['Covariance'])
1596            line = '  ' + PutInCol(G2mth.ValEsd(T,-1.),6)
1597            for i in terms:
1598                line += PutInCol(G2mth.ValEsd(cellList[i],cellSig[i]),12)
1599            line += ' ' + datablockidDict[h]
1600            WriteCIFitem(self.fp, line)
1601        for h,T in DijTlist.items():
1602            pId = phasedict['pId']
1603            hId = self.Histograms[h]['hId']
1604            cellList,cellSig = G2stIO.getCellSU(pId,hId,
1605                                        phasedict['General']['SGData'],
1606                                        self.parmDict,
1607                                        self.OverallParms['Covariance'])
1608            line = '  ' + PutInCol(G2mth.ValEsd(T,-1.),6)
1609            for i in terms:
1610                line += PutInCol(G2mth.ValEsd(cellList[i],cellSig[i]),12)
1611            line += ' ' + datablockidDict[h]
1612            WriteCIFitem(self.fp, line)       
1613
1614    def _Exporter(self,event=None,phaseOnly=None,histOnly=None):
1615        '''Basic code to export a CIF. Export can be full or simple, as set by
1616        phaseOnly and histOnly which skips distances & angles, etc.
1617
1618        :param bool phaseOnly: used to export only one phase
1619        :param bool histOnly: used to export only one histogram
1620        '''
1621
1622#***** define functions for export method =======================================
1623        def WriteAudit():
1624            'Write the CIF audit values. Perhaps should be in a single element loop.'
1625            WriteCIFitem(self.fp, '_audit_creation_method',
1626                         'created in GSAS-II')
1627            WriteCIFitem(self.fp, '_audit_creation_date',self.CIFdate)
1628            if self.author:
1629                WriteCIFitem(self.fp, '_audit_author_name',self.author)
1630            WriteCIFitem(self.fp, '_audit_update_record',
1631                         self.CIFdate+'  Initial software-generated CIF')
1632
1633        def WriteOverall(mode=None):
1634            '''Write out overall refinement information.
1635
1636            More could be done here, but this is a good start.
1637            '''
1638            if self.ifPWDR:
1639                WriteCIFitem(self.fp, '_pd_proc_info_datetime', self.CIFdate)
1640                WriteCIFitem(self.fp, '_pd_calc_method', 'Rietveld Refinement')
1641               
1642            #WriteCIFitem(self.fp, '_refine_ls_shift/su_mean',DAT2)
1643            WriteCIFitem(self.fp, '_computing_structure_refinement','GSAS-II (Toby & Von Dreele, J. Appl. Cryst. 46, 544-549, 2013)')
1644            if self.ifHKLF:
1645                controls = self.OverallParms['Controls']
1646                try:
1647                    if controls['F**2']:
1648                        thresh = 'F**2>%.1fu(F**2)'%(controls['UsrReject']['minF/sig'])
1649                    else:
1650                        thresh = 'F>%.1fu(F)'%(controls['UsrReject']['minF/sig'])
1651                    WriteCIFitem(self.fp, '_reflns_threshold_expression', thresh)
1652                except KeyError:
1653                    pass
1654            WriteCIFitem(self.fp, '_refine_ls_matrix_type','full')
1655
1656            if mode == 'seq': return
1657            try:
1658                vars = str(len(self.OverallParms['Covariance']['varyList']))
1659            except:
1660                vars = '?'
1661            WriteCIFitem(self.fp, '_refine_ls_number_parameters',vars)
1662            try:
1663                GOF = G2mth.ValEsd(self.OverallParms['Covariance']['Rvals']['GOF'],-0.009)
1664            except:
1665                GOF = '?'
1666            WriteCIFitem(self.fp, '_refine_ls_goodness_of_fit_all',GOF)
1667            DAT1 = self.OverallParms['Covariance']['Rvals'].get('Max shft/sig',0.0)
1668            if DAT1:
1669                WriteCIFitem(self.fp, '_refine_ls_shift/su_max','%.4f'%DAT1)
1670
1671            # get restraint info
1672            # restraintDict = self.OverallParms.get('Restraints',{})
1673            # for i in  self.OverallParms['Constraints']:
1674            #     print i
1675            #     for j in self.OverallParms['Constraints'][i]:
1676            #         print j
1677            #WriteCIFitem(self.fp, '_refine_ls_number_restraints',TEXT)
1678            # other things to consider reporting
1679            # _refine_ls_number_reflns
1680            # _refine_ls_goodness_of_fit_obs
1681            # _refine_ls_wR_factor_obs
1682            # _refine_ls_restrained_S_all
1683            # _refine_ls_restrained_S_obs
1684
1685            # include an overall profile r-factor, if there is more than one powder histogram
1686            R = '%.5f'%(self.OverallParms['Covariance']['Rvals']['Rwp']/100.)
1687            WriteCIFitem(self.fp, '\n# OVERALL WEIGHTED R-FACTOR')
1688            WriteCIFitem(self.fp, '_refine_ls_wR_factor_obs',R)
1689                # _refine_ls_R_factor_all
1690                # _refine_ls_R_factor_obs
1691            #WriteCIFitem(self.fp, '_refine_ls_matrix_type','userblocks')
1692
1693        def WriteOverallMM(mode=None):
1694            '''Write out overall refinement information.
1695
1696            More could be done here, but this is a good start.
1697            '''
1698            if self.ifPWDR:
1699                WriteCIFitem(self.fp, '_pd_proc_info_datetime', self.CIFdate)
1700                WriteCIFitem(self.fp, '_pd_calc_method', 'Rietveld Refinement')
1701               
1702            #WriteCIFitem(self.fp, '_refine.ls_shift_over_su_mean',DAT2)
1703            WriteCIFitem(self.fp, '_computing_structure_refinement','GSAS-II (Toby & Von Dreele, J. Appl. Cryst. 46, 544-549, 2013)')
1704            if self.ifHKLF:
1705                controls = self.OverallParms['Controls']
1706                try:
1707                    if controls['F**2']:
1708                        thresh = 'F**2>%.1fu(F**2)'%(controls['UsrReject']['minF/sig'])
1709                    else:
1710                        thresh = 'F>%.1fu(F)'%(controls['UsrReject']['minF/sig'])
1711                    WriteCIFitem(self.fp, '_reflns.threshold_expression', thresh)
1712                except KeyError:
1713                    pass
1714            WriteCIFitem(self.fp, '_refine.ls_matrix_type','full')
1715
1716            if mode == 'seq': return
1717            try:
1718                vars = str(len(self.OverallParms['Covariance']['varyList']))
1719            except:
1720                vars = '?'
1721            WriteCIFitem(self.fp, '_refine.ls_number_parameters',vars)
1722            try:
1723                GOF = G2mth.ValEsd(self.OverallParms['Covariance']['Rvals']['GOF'],-0.009)
1724            except:
1725                GOF = '?'
1726            WriteCIFitem(self.fp, '_refine.ls_goodness_of_fit_all',GOF)
1727            DAT1 = self.OverallParms['Covariance']['Rvals'].get('Max shft/sig',0.0)
1728            if DAT1:
1729                WriteCIFitem(self.fp, '_refine.ls_shift_over_su_max','%.4f'%DAT1)
1730
1731            # get restraint info
1732            # restraintDict = self.OverallParms.get('Restraints',{})
1733            # for i in  self.OverallParms['Constraints']:
1734            #     print i
1735            #     for j in self.OverallParms['Constraints'][i]:
1736            #         print j
1737            #WriteCIFitem(self.fp, '_refine_ls_number_restraints',TEXT)
1738            # other things to consider reporting
1739            # _refine_ls_number_reflns
1740            # _refine_ls_goodness_of_fit_obs
1741            # _refine_ls_wR_factor_obs
1742            # _refine_ls_restrained_S_all
1743            # _refine_ls_restrained_S_obs
1744
1745            # include an overall profile r-factor, if there is more than one powder histogram
1746            R = '%.5f'%(self.OverallParms['Covariance']['Rvals']['Rwp']/100.)
1747            WriteCIFitem(self.fp, '\n# OVERALL WEIGHTED R-FACTOR')
1748            WriteCIFitem(self.fp, '_refine.ls_wR_factor_obs',R)
1749
1750        def writeCIFtemplate(G2dict,tmplate,defaultname='',
1751                                 cifKey="CIF_template"):
1752            '''Write out the selected or edited CIF template
1753            An unedited CIF template file is copied, comments intact; an edited
1754            CIF template is written out from PyCifRW which of course strips comments.
1755            In all cases the initial data_ header is stripped (there should only be one!)
1756            '''
1757            CIFobj = G2dict.get(cifKey)
1758            if defaultname:
1759                defaultname = G2obj.StripUnicode(defaultname)
1760                defaultname = re.sub(r'[^a-zA-Z0-9_-]','',defaultname)
1761                defaultname = tmplate + "_" + defaultname + ".cif"
1762            else:
1763                defaultname = ''
1764            templateDefName = 'template_'+tmplate+'.cif'
1765            if not CIFobj: # copying a template
1766                lbl = 'Standard version'
1767                for pth in [os.getcwd()]+sys.path:
1768                    fil = os.path.join(pth,defaultname)
1769                    if os.path.exists(fil) and defaultname: break
1770                else:
1771                    for pth in sys.path:
1772                        fil = os.path.join(pth,templateDefName)
1773                        if os.path.exists(fil): break
1774                    else:
1775                        print(CIFobj+' not found in path!')
1776                        return
1777                fp = open(fil,'r')
1778                txt = fp.read()
1779                fp.close()
1780            elif type(CIFobj) is not list and type(CIFobj) is not tuple:
1781                lbl = 'Saved version'
1782                if not os.path.exists(CIFobj):
1783                    print("Error: requested template file has disappeared: "+CIFobj)
1784                    return
1785                fp = open(CIFobj,'r')
1786                txt = fp.read()
1787                fp.close()
1788            else:
1789                lbl = 'Project-specific version'
1790                txt = dict2CIF(CIFobj[0],CIFobj[1]).WriteOut()
1791            # remove the PyCifRW header, if present
1792            #if txt.find('PyCifRW') > -1 and txt.find('data_') > -1:
1793            pre = txt.index("data_")
1794            restofline = txt.index("\n",pre)
1795            name = txt[pre+5:restofline]
1796            txt = "\n# {} of {} template follows{}".format(
1797                lbl, name, txt[restofline:])
1798            #txt = txt.replace('data_','#')
1799            WriteCIFitem(self.fp, txt)
1800
1801        def FormatSH(phasenam):
1802            'Format a full spherical harmonics texture description as a string'
1803            phasedict = self.Phases[phasenam] # pointer to current phase info
1804            pfx = str(phasedict['pId'])+'::'
1805            s = ""
1806            textureData = phasedict['General']['SH Texture']
1807            if textureData.get('Order'):
1808                s += "Spherical Harmonics correction. Order = "+str(textureData['Order'])
1809                s += " Model: " + str(textureData['Model']) + "\n    Orientation angles: "
1810                for name in ['omega','chi','phi']:
1811                    aname = pfx+'SH '+name
1812                    s += name + " = "
1813                    sig = self.sigDict.get(aname,-0.09)
1814                    s += G2mth.ValEsd(self.parmDict[aname],sig)
1815                    s += "; "
1816                s += "\n"
1817                s1 = "    Coefficients:  "
1818                for name in textureData['SH Coeff'][1]:
1819                    aname = pfx+name
1820                    if len(s1) > 60:
1821                        s += s1 + "\n"
1822                        s1 = "    "
1823                    s1 += aname + ' = '
1824                    sig = self.sigDict.get(aname,-0.0009)
1825                    s1 += G2mth.ValEsd(self.parmDict[aname],sig)
1826                    s1 += "; "
1827                s += s1
1828            return s
1829
1830        def FormatHAPpo(phasenam):
1831            '''return the March-Dollase/SH correction for every
1832            histogram in the current phase formatted into a
1833            character string
1834            '''
1835            phasedict = self.Phases[phasenam] # pointer to current phase info
1836            s = ''
1837            for histogram in sorted(phasedict['Histograms']):
1838                if histogram.startswith("HKLF"): continue # powder only
1839                if not self.Phases[phasenam]['Histograms'][histogram]['Use']: continue
1840                Histogram = self.Histograms.get(histogram)
1841                if not Histogram: continue
1842                hapData = phasedict['Histograms'][histogram]
1843                if hapData['Pref.Ori.'][0] == 'MD':
1844                    aname = str(phasedict['pId'])+':'+str(Histogram['hId'])+':MD'
1845                    if self.parmDict.get(aname,1.0) != 1.0: continue
1846                    sig = self.sigDict.get(aname,-0.009)
1847                    if s != "": s += '\n'
1848                    s += 'March-Dollase correction'
1849                    if len(self.powderDict) > 1:
1850                        s += ', histogram '+str(Histogram['hId']+1)
1851                    s += ' coef. = ' + G2mth.ValEsd(self.parmDict[aname],sig)
1852                    s += ' axis = ' + str(hapData['Pref.Ori.'][3])
1853                else: # must be SH
1854                    if s != "": s += '\n'
1855                    s += 'Simple spherical harmonic correction'
1856                    if len(self.powderDict) > 1:
1857                        s += ', histogram '+str(Histogram['hId']+1)
1858                    s += ' Order = '+str(hapData['Pref.Ori.'][4])+'\n'
1859                    s1 = "    Coefficients:  "
1860                    for item in hapData['Pref.Ori.'][5]:
1861                        aname = str(phasedict['pId'])+':'+str(Histogram['hId'])+':'+item
1862                        if len(s1) > 60:
1863                            s += s1 + "\n"
1864                            s1 = "    "
1865                        s1 += aname + ' = '
1866                        sig = self.sigDict.get(aname,-0.0009)
1867                        s1 += G2mth.ValEsd(self.parmDict[aname],sig)
1868                        s1 += "; "
1869                    s += s1
1870            return s
1871
1872        def FormatBackground(bkg,hId):
1873            '''Display the Background information as a descriptive text string.
1874
1875            TODO: this needs to be expanded to show the diffuse peak and
1876            Debye term information as well. (Bob)
1877
1878            :returns: the text description (str)
1879            '''
1880            hfx = ':'+str(hId)+':'
1881            fxn, bkgdict = bkg
1882            terms = fxn[2]
1883            txt = 'Background function: "'+fxn[0]+'" function with '+str(terms)+' terms:\n'
1884            l = "    "
1885            for i,v in enumerate(fxn[3:]):
1886                name = '%sBack;%d'%(hfx,i)
1887                sig = self.sigDict.get(name,-0.009)
1888                if len(l) > 60:
1889                    txt += l + '\n'
1890                    l = '    '
1891                l += G2mth.ValEsd(v,sig)+', '
1892            txt += l
1893            if bkgdict['nDebye']:
1894                txt += '\n  Background Debye function parameters: A, R, U:'
1895                names = ['A;','R;','U;']
1896                for i in range(bkgdict['nDebye']):
1897                    txt += '\n    '
1898                    for j in range(3):
1899                        name = hfx+'Debye'+names[j]+str(i)
1900                        sig = self.sigDict.get(name,-0.009)
1901                        txt += G2mth.ValEsd(bkgdict['debyeTerms'][i][2*j],sig)+', '
1902            if bkgdict['nPeaks']:
1903                txt += '\n  Background peak parameters: pos, int, sig, gam:'
1904                names = ['pos;','int;','sig;','gam;']
1905                for i in range(bkgdict['nPeaks']):
1906                    txt += '\n    '
1907                    for j in range(4):
1908                        name = hfx+'BkPk'+names[j]+str(i)
1909                        sig = self.sigDict.get(name,-0.009)
1910                        txt += G2mth.ValEsd(bkgdict['peaksList'][i][2*j],sig)+', '
1911            return txt
1912
1913        def FormatInstProfile(instparmdict,hId):
1914            '''Format the instrumental profile parameters with a
1915            string description. Will only be called on PWDR histograms
1916            '''
1917            s = ''
1918            inst = instparmdict[0]
1919            hfx = ':'+str(hId)+':'
1920            if 'C' in inst['Type'][0]:
1921                s = 'Finger-Cox-Jephcoat function parameters U, V, W, X, Y, SH/L:\n'
1922                s += '  peak variance(Gauss) = Utan(Th)^2^+Vtan(Th)+W:\n'
1923                s += '  peak HW(Lorentz) = X/cos(Th)+Ytan(Th); SH/L = S/L+H/L\n'
1924                s += '  U, V, W in (centideg)^2^, X & Y in centideg\n    '
1925                for item in ['U','V','W','X','Y','SH/L']:
1926                    name = hfx+item
1927                    sig = self.sigDict.get(name,-0.009)
1928                    s += G2mth.ValEsd(inst[item][1],sig)+', '
1929            elif 'T' in inst['Type'][0]:    #to be tested after TOF Rietveld done
1930                s = 'Von Dreele-Jorgenson-Windsor function parameters\n'+ \
1931                    '   alpha, beta-0, beta-1, beta-q, sig-0, sig-1, sig-2, sig-q, X, Y:\n    '
1932                for item in ['alpha','beta-0','beta-1','beta-q','sig-0','sig-1','sig-2','sig-q','X','Y']:
1933                    name = hfx+item
1934                    sig = self.sigDict.get(name,-0.009)
1935                    s += G2mth.ValEsd(inst[item][1],sig)+', '
1936            return s
1937
1938        def FormatPhaseProfile(phasenam,hist=''):
1939            '''Format the phase-related profile parameters (size/strain)
1940            with a string description.
1941            return an empty string or None if there are no
1942            powder histograms for this phase.
1943            '''
1944            s = ''
1945            phasedict = self.Phases[phasenam] # pointer to current phase info
1946            if hist:
1947                parmDict = self.seqData[hist]['parmDict']
1948                sigDict = dict(zip(self.seqData[hist]['varyList'],self.seqData[hist]['sig']))
1949            else:
1950                parmDict = self.parmDict
1951                sigDict = self.sigDict
1952           
1953            SGData = phasedict['General'] ['SGData']
1954            for histogram in sorted(phasedict['Histograms']):
1955                if hist is not None and hist != histogram: continue
1956                if histogram.startswith("HKLF"): continue # powder only
1957                Histogram = self.Histograms.get(histogram)
1958                if not Histogram: continue
1959                hapData = phasedict['Histograms'][histogram]
1960                pId = phasedict['pId']
1961                hId = Histogram['hId']
1962                phfx = '%d:%d:'%(pId,hId)
1963                size = hapData['Size']
1964                mustrain = hapData['Mustrain']
1965                hstrain = hapData['HStrain']
1966                if s: s += '\n'
1967                if len(self.powderDict) > 1: # if one histogram, no ambiguity
1968                    s += '  Parameters for histogram #{:} {:} & phase {:}\n'.format(
1969                        str(hId),str(histogram),phasenam)
1970                s += '  Crystallite size in microns with "%s" model:\n  '%(size[0])
1971                names = ['Size;i','Size;mx']
1972                if 'uniax' in size[0]:
1973                    names = ['Size;i','Size;a','Size;mx']
1974                    s += 'anisotropic axis is %s\n  '%(str(size[3]))
1975                    s += 'parameters: equatorial size, axial size, G/L mix\n    '
1976                    for i,item in enumerate(names):
1977                        name = phfx+item
1978                        val = parmDict.get(name,size[1][i])
1979                        sig = sigDict.get(name,-0.009)
1980                        s += G2mth.ValEsd(val,sig)+', '
1981                elif 'ellip' in size[0]:
1982                    s += 'parameters: S11, S22, S33, S12, S13, S23, G/L mix\n    '
1983                    for i in range(6):
1984                        name = phfx+'Size:'+str(i)
1985                        val = parmDict.get(name,size[4][i])
1986                        sig = sigDict.get(name,-0.009)
1987                        s += G2mth.ValEsd(val,sig)+', '
1988                    sig = sigDict.get(phfx+'Size;mx',-0.009)
1989                    s += G2mth.ValEsd(size[1][2],sig)+', '
1990                else:       #isotropic
1991                    s += 'parameters: Size, G/L mix\n    '
1992                    i = 0
1993                    for item in names:
1994                        name = phfx+item
1995                        val = parmDict.get(name,size[1][i])
1996                        sig = sigDict.get(name,-0.009)
1997                        s += G2mth.ValEsd(val,sig)+', '
1998                        i = 2    #skip the aniso value
1999                s += '\n  Microstrain, "%s" model (10^6^ * delta Q/Q)\n  '%(mustrain[0])
2000                names = ['Mustrain;i','Mustrain;mx']
2001                if 'uniax' in mustrain[0]:
2002                    names = ['Mustrain;i','Mustrain;a','Mustrain;mx']
2003                    s += 'anisotropic axis is %s\n  '%(str(size[3]))
2004                    s += 'parameters: equatorial mustrain, axial mustrain, G/L mix\n    '
2005                    for i,item in enumerate(names):
2006                        name = phfx+item
2007                        val = parmDict.get(name,mustrain[1][i])
2008                        sig = sigDict.get(name,-0.009)
2009                        s += G2mth.ValEsd(val,sig)+', '
2010                elif 'general' in mustrain[0]:
2011                    names = 'parameters: '
2012                    for i,name in enumerate(G2spc.MustrainNames(SGData)):
2013                        names += name+', '
2014                        if i == 9:
2015                            names += '\n  '
2016                    names += 'G/L mix\n    '
2017                    s += names
2018                    txt = ''
2019                    for i in range(len(mustrain[4])):
2020                        name = phfx+'Mustrain:'+str(i)
2021                        val = parmDict.get(name,mustrain[4][i])
2022                        sig = sigDict.get(name,-0.009)
2023                        if len(txt) > 60:
2024                            s += txt+'\n    '
2025                            txt = ''
2026                        txt += G2mth.ValEsd(val,sig)+', '
2027                    s += txt
2028                    name = phfx+'Mustrain;mx'
2029                    val = parmDict.get(name,mustrain[1][2])
2030                    sig = sigDict.get(name,-0.009)
2031                    s += G2mth.ValEsd(val,sig)+', '
2032
2033                else:       #isotropic
2034                    s += '  parameters: Mustrain, G/L mix\n    '
2035                    i = 0
2036                    for item in names:
2037                        name = phfx+item
2038                        val = parmDict.get(name,mustrain[1][i])
2039                        sig = sigDict.get(name,-0.009)
2040                        s += G2mth.ValEsd(val,sig)+', '
2041                        i = 2    #skip the aniso value
2042                s1 = \n  Macrostrain parameters: '
2043                names = G2spc.HStrainNames(SGData)
2044                for name in names:
2045                    s1 += name+', '
2046                s1 += '\n    '
2047                macrostrain = False
2048                for i in range(len(names)):
2049                    name = phfx+names[i]
2050                    val = parmDict.get(name,hstrain[0][i])
2051                    sig = sigDict.get(name,-0.000009)
2052                    s1 += G2mth.ValEsd(val,sig)+', '
2053                    if hstrain[0][i]: macrostrain = True
2054                if macrostrain:
2055                    s += s1 + '\n'
2056                    # show revised lattice parameters here someday
2057                else:
2058                    s += '\n'
2059            return s
2060
2061        def MakeUniqueLabel(lbl,labellist):
2062            'Make sure that every atom label is unique'
2063            lbl = lbl.strip()
2064            if not lbl: # deal with a blank label
2065                lbl = 'A_1'
2066            if lbl not in labellist:
2067                labellist.append(lbl)
2068                return lbl
2069            i = 1
2070            prefix = lbl
2071            if '_' in lbl:
2072                prefix = lbl[:lbl.rfind('_')]
2073                suffix = lbl[lbl.rfind('_')+1:]
2074                try:
2075                    i = int(suffix)+1
2076                except:
2077                    pass
2078            while prefix+'_'+str(i) in labellist:
2079                i += 1
2080            else:
2081                lbl = prefix+'_'+str(i)
2082                labellist.append(lbl)
2083
2084        def WriteDistances(phasenam):
2085            '''Report bond distances and angles for the CIF
2086
2087            Note that _geom_*_symmetry_* fields are values of form
2088            n_klm where n is the symmetry operation in SymOpList (counted
2089            starting with 1) and (k-5, l-5, m-5) are translations to add
2090            to (x,y,z). See
2091            http://www.iucr.org/__data/iucr/cifdic_html/1/cif_core.dic/Igeom_angle_site_symmetry_.html
2092
2093            TODO: need a method to select publication flags for distances/angles
2094            '''
2095            phasedict = self.Phases[phasenam] # pointer to current phase info
2096            Atoms = phasedict['Atoms']
2097            generalData = phasedict['General']
2098            # create a dict for storing Pub flag for bonds/angles, if needed
2099            if phasedict['General'].get("DisAglHideFlag") is None:
2100                phasedict['General']["DisAglHideFlag"] = {}
2101            DisAngSel = phasedict['General']["DisAglHideFlag"]
2102            cx,ct,cs,cia = phasedict['General']['AtomPtrs']
2103            cn = ct-1
2104            fpfx = str(phasedict['pId'])+'::Afrac:'
2105            cfrac = cx+3
2106            DisAglData = {}
2107            # create a list of atoms, but skip atoms with zero occupancy
2108            xyz = []
2109            fpfx = str(phasedict['pId'])+'::Afrac:'
2110            for i,atom in enumerate(Atoms):
2111                if self.parmDict.get(fpfx+str(i),atom[cfrac]) == 0.0: continue
2112                xyz.append([i,]+atom[cn:cn+2]+atom[cx:cx+3])
2113            if 'DisAglCtls' not in generalData:
2114                # should not happen, since DisAglDialog should be called
2115                # for all phases before getting here
2116                dlg = G2G.DisAglDialog(
2117                    self.G2frame,
2118                    {},
2119                    generalData)
2120                if dlg.ShowModal() == wx.ID_OK:
2121                    generalData['DisAglCtls'] = dlg.GetData()
2122                else:
2123                    dlg.Destroy()
2124                    return
2125                dlg.Destroy()
2126            DisAglData['OrigAtoms'] = xyz
2127            DisAglData['TargAtoms'] = xyz
2128            SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
2129                generalData['SGData'])
2130
2131#            xpandSGdata = generalData['SGData'].copy()
2132#            xpandSGdata.update({'SGOps':symOpList,
2133#                                'SGInv':False,
2134#                                'SGLatt':'P',
2135#                                'SGCen':np.array([[0, 0, 0]]),})
2136#            DisAglData['SGData'] = xpandSGdata
2137            DisAglData['SGData'] = generalData['SGData'].copy()
2138
2139            DisAglData['Cell'] = generalData['Cell'][1:] #+ volume
2140            if 'pId' in phasedict:
2141                DisAglData['pId'] = phasedict['pId']
2142                DisAglData['covData'] = self.OverallParms['Covariance']
2143            try:
2144                AtomLabels,DistArray,AngArray = G2stMn.RetDistAngle(
2145                    generalData['DisAglCtls'],
2146                    DisAglData)
2147            except KeyError:        # inside DistAngle for missing atom types in DisAglCtls
2148                print(u'**** ERROR computing distances & angles for phase {} ****\nresetting to default values'.format(phasenam))
2149                data = generalData['DisAglCtls'] = {}
2150                data['Name'] = generalData['Name']
2151                data['Factors'] = [0.85,0.85]
2152                data['AtomTypes'] = generalData['AtomTypes']
2153                data['BondRadii'] = generalData['BondRadii'][:]
2154                data['AngleRadii'] = generalData['AngleRadii'][:]
2155                try:
2156                    AtomLabels,DistArray,AngArray = G2stMn.RetDistAngle(
2157                        generalData['DisAglCtls'],
2158                        DisAglData)
2159                except:
2160                    print('Reset failed. To fix this, use the Reset button in the "edit distance/angle menu" for this phase')
2161                    return
2162
2163            # loop over interatomic distances for this phase
2164            WriteCIFitem(self.fp, '\n# MOLECULAR GEOMETRY')
2165            First = True
2166            for i in sorted(AtomLabels.keys()):
2167                Dist = DistArray[i]
2168                for D in Dist:
2169                    line = '  '+PutInCol(AtomLabels[i],6)+PutInCol(AtomLabels[D[0]],6)
2170                    sig = D[4]
2171                    if sig == 0: sig = -0.00009
2172                    line += PutInCol(G2mth.ValEsd(D[3],sig,True),10)
2173                    line += "  1_555 "
2174                    symopNum = G2opcodes.index(D[2])
2175                    line += " {:3d}_".format(symopNum+1)
2176                    for d,o in zip(D[1],offsetList[symopNum]):
2177                        line += "{:1d}".format(d-o+5)
2178                    if DisAngSel.get((i,tuple(D[0:3]))):
2179                        line += " no"
2180                    else:
2181                        line += " yes"
2182                    if First:
2183                        First = False
2184                        WriteCIFitem(self.fp, 'loop_' +
2185                         '\n   _geom_bond_atom_site_label_1' +
2186                         '\n   _geom_bond_atom_site_label_2' +
2187                         '\n   _geom_bond_distance' +
2188                         '\n   _geom_bond_site_symmetry_1' +
2189                         '\n   _geom_bond_site_symmetry_2' +
2190                         '\n   _geom_bond_publ_flag')
2191                    WriteCIFitem(self.fp, line)
2192
2193            # loop over interatomic angles for this phase
2194            First = True
2195            for i in sorted(AtomLabels.keys()):
2196                Dist = DistArray[i]
2197                for k,j,tup in AngArray[i]:
2198                    Dj = Dist[j]
2199                    Dk = Dist[k]
2200                    line = '  '+PutInCol(AtomLabels[Dj[0]],6)+PutInCol(AtomLabels[i],6)+PutInCol(AtomLabels[Dk[0]],6)
2201                    sig = tup[1]
2202                    if sig == 0: sig = -0.009
2203                    line += PutInCol(G2mth.ValEsd(tup[0],sig,True),10)
2204                    line += " {:3d}_".format(G2opcodes.index(Dj[2])+1)
2205                    for d in Dj[1]:
2206                        line += "{:1d}".format(d+5)
2207                    line += "  1_555 "
2208                    line += " {:3d}_".format(G2opcodes.index(Dk[2])+1)
2209                    for d in Dk[1]:
2210                        line += "{:1d}".format(d+5)
2211                    key = (tuple(Dk[0:3]),i,tuple(Dj[0:3]))
2212                    if DisAngSel.get(key):
2213                        line += " no"
2214                    else:
2215                        line += " yes"
2216                    if First:
2217                        First = False
2218                        WriteCIFitem(self.fp, '\nloop_' +
2219                         '\n   _geom_angle_atom_site_label_1' +
2220                         '\n   _geom_angle_atom_site_label_2' +
2221                         '\n   _geom_angle_atom_site_label_3' +
2222                         '\n   _geom_angle' +
2223                         '\n   _geom_angle_site_symmetry_1' +
2224                         '\n   _geom_angle_site_symmetry_2' +
2225                         '\n   _geom_angle_site_symmetry_3' +
2226                         '\n   _geom_angle_publ_flag')
2227                    WriteCIFitem(self.fp, line)
2228
2229
2230        def WriteSeqDistances(phasenam,histname,phasedict,cellList,seqData):
2231            '''Report bond distances and angles for the CIF from a Sequential fit
2232
2233            Note that _geom_*_symmetry_* fields are values of form
2234            n_klm where n is the symmetry operation in SymOpList (counted
2235            starting with 1) and (k-5, l-5, m-5) are translations to add
2236            to (x,y,z). See
2237            http://www.iucr.org/__data/iucr/cifdic_html/1/cif_core.dic/Igeom_angle_site_symmetry_.html
2238
2239            TODO: this is based on WriteDistances and could likely be merged with that
2240            without too much work. Note also that G2stMn.RetDistAngle is pretty slow for
2241            sequential fits, since it is called so many times.
2242            '''
2243            #breakpoint()
2244            Atoms = phasedict['Atoms']
2245            generalData = phasedict['General']
2246            parmDict = seqData[histname]['parmDict']
2247#            sigDict = dict(zip(seqData[hist]['varyList'],seqData[hist]['sig']))
2248            # create a dict for storing Pub flag for bonds/angles, if needed
2249            if phasedict['General'].get("DisAglHideFlag") is None:
2250                phasedict['General']["DisAglHideFlag"] = {}
2251            DisAngSel = phasedict['General']["DisAglHideFlag"]
2252            cx,ct,cs,cia = phasedict['General']['AtomPtrs']
2253            cn = ct-1
2254#            fpfx = str(phasedict['pId'])+'::Afrac:'
2255            cfrac = cx+3
2256            DisAglData = {}
2257            # create a list of atoms, but skip atoms with zero occupancy
2258            xyz = []
2259            fpfx = str(phasedict['pId'])+'::Afrac:'
2260            for i,atom in enumerate(Atoms):
2261                if parmDict.get(fpfx+str(i),atom[cfrac]) == 0.0: continue
2262                thisatom = [i] + atom[cn:cn+2]
2263                for j,lab in enumerate(['x','y','z']):
2264                    xyzkey = str(phasedict['pId'])+'::A'+ lab + ':' +str(i)
2265                    thisatom.append(parmDict.get(xyzkey,atom[cx+j]))
2266                xyz.append(thisatom)
2267            DisAglData['OrigAtoms'] = xyz
2268            DisAglData['TargAtoms'] = xyz
2269            SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
2270                generalData['SGData'])
2271
2272#            xpandSGdata = generalData['SGData'].copy()
2273#            xpandSGdata.update({'SGOps':symOpList,
2274#                                'SGInv':False,
2275#                                'SGLatt':'P',
2276#                                'SGCen':np.array([[0, 0, 0]]),})
2277#            DisAglData['SGData'] = xpandSGdata
2278            DisAglData['SGData'] = generalData['SGData'].copy()
2279
2280            DisAglData['Cell'] = cellList  #+ volume
2281            if 'pId' in phasedict:
2282                DisAglData['pId'] = phasedict['pId']
2283                DisAglData['covData'] = seqData[histname]
2284                # self.OverallParms['Covariance']
2285            try:
2286                AtomLabels,DistArray,AngArray = G2stMn.RetDistAngle(
2287                    generalData['DisAglCtls'],
2288                    DisAglData)
2289            except KeyError:        # inside DistAngle for missing atom types in DisAglCtls
2290                print(u'**** ERROR computing distances & angles for phase {} ****\nresetting to default values'.format(phasenam))
2291                data = generalData['DisAglCtls'] = {}
2292                data['Name'] = generalData['Name']
2293                data['Factors'] = [0.85,0.85]
2294                data['AtomTypes'] = generalData['AtomTypes']
2295                data['BondRadii'] = generalData['BondRadii'][:]
2296                data['AngleRadii'] = generalData['AngleRadii'][:]
2297                try:
2298                    AtomLabels,DistArray,AngArray = G2stMn.RetDistAngle(
2299                        generalData['DisAglCtls'],
2300                        DisAglData)
2301                except:
2302                    print('Reset failed. To fix this, use the Reset button in the "edit distance/angle menu" for this phase')
2303                    return
2304
2305            # loop over interatomic distances for this phase
2306            WriteCIFitem(self.fp, '\n# MOLECULAR GEOMETRY')
2307            First = True
2308            for i in sorted(AtomLabels.keys()):
2309                Dist = DistArray[i]
2310                for D in Dist:
2311                    line = '  '+PutInCol(AtomLabels[i],6)+PutInCol(AtomLabels[D[0]],6)
2312                    sig = D[4]
2313                    if sig == 0: sig = -0.00009
2314                    line += PutInCol(G2mth.ValEsd(D[3],sig,True),10)
2315                    line += "  1_555 "
2316                    symopNum = G2opcodes.index(D[2])
2317                    line += " {:3d}_".format(symopNum+1)
2318                    for d,o in zip(D[1],offsetList[symopNum]):
2319                        line += "{:1d}".format(d-o+5)
2320                    if DisAngSel.get((i,tuple(D[0:3]))):
2321                        line += " no"
2322                    else:
2323                        line += " yes"
2324                    if First:
2325                        First = False
2326                        WriteCIFitem(self.fp, 'loop_' +
2327                         '\n   _geom_bond_atom_site_label_1' +
2328                         '\n   _geom_bond_atom_site_label_2' +
2329                         '\n   _geom_bond_distance' +
2330                         '\n   _geom_bond_site_symmetry_1' +
2331                         '\n   _geom_bond_site_symmetry_2' +
2332                         '\n   _geom_bond_publ_flag')
2333                    WriteCIFitem(self.fp, line)
2334
2335            # loop over interatomic angles for this phase
2336            First = True
2337            for i in sorted(AtomLabels.keys()):
2338                Dist = DistArray[i]
2339                for k,j,tup in AngArray[i]:
2340                    Dj = Dist[j]
2341                    Dk = Dist[k]
2342                    line = '  '+PutInCol(AtomLabels[Dj[0]],6)+PutInCol(AtomLabels[i],6)+PutInCol(AtomLabels[Dk[0]],6)
2343                    sig = tup[1]
2344                    if sig == 0: sig = -0.009
2345                    line += PutInCol(G2mth.ValEsd(tup[0],sig,True),10)
2346                    line += " {:3d}_".format(G2opcodes.index(Dj[2])+1)
2347                    for d in Dj[1]:
2348                        line += "{:1d}".format(d+5)
2349                    line += "  1_555 "
2350                    line += " {:3d}_".format(G2opcodes.index(Dk[2])+1)
2351                    for d in Dk[1]:
2352                        line += "{:1d}".format(d+5)
2353                    key = (tuple(Dk[0:3]),i,tuple(Dj[0:3]))
2354                    if DisAngSel.get(key):
2355                        line += " no"
2356                    else:
2357                        line += " yes"
2358                    if First:
2359                        First = False
2360                        WriteCIFitem(self.fp, '\nloop_' +
2361                         '\n   _geom_angle_atom_site_label_1' +
2362                         '\n   _geom_angle_atom_site_label_2' +
2363                         '\n   _geom_angle_atom_site_label_3' +
2364                         '\n   _geom_angle' +
2365                         '\n   _geom_angle_site_symmetry_1' +
2366                         '\n   _geom_angle_site_symmetry_2' +
2367                         '\n   _geom_angle_site_symmetry_3' +
2368                         '\n   _geom_angle_publ_flag')
2369                    WriteCIFitem(self.fp, line)
2370
2371        def WriteSeqOverallPhaseInfo(phasenam,histblk):
2372            'Write out the phase information for the selected phase for the overall block in a sequential fit'
2373            WriteCIFitem(self.fp, '# overall phase info for '+str(phasenam) + ' follows')
2374            phasedict = self.Phases[phasenam] # pointer to current phase info
2375            WriteCIFitem(self.fp, '_pd_phase_name', phasenam)
2376
2377            WriteCIFitem(self.fp, '_symmetry_cell_setting',
2378                         phasedict['General']['SGData']['SGSys'])
2379
2380              # moved to WriteSeqPhaseVals()
2381#             if phasedict['General']['Type'] in ['nuclear','macromolecular']:
2382#                 spacegroup = phasedict['General']['SGData']['SpGrp'].strip()
2383#                 # regularize capitalization and remove trailing H/R
2384#                 spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
2385#                 WriteCIFitem(self.fp, '_symmetry_space_group_name_H-M',spacegroup)
2386   
2387#                 # generate symmetry operations including centering and center of symmetry
2388#                 SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
2389#                     phasedict['General']['SGData'])
2390#                 WriteCIFitem(self.fp, 'loop_\n    _space_group_symop_id\n    _space_group_symop_operation_xyz')
2391#                 for i,op in enumerate(SymOpList,start=1):
2392#                     WriteCIFitem(self.fp, '   {:3d}  {:}'.format(i,op.lower()))
2393#             elif phasedict['General']['Type'] == 'magnetic':
2394#                 parentSpGrp = phasedict['General']['SGData']['SpGrp'].strip()
2395#                 parentSpGrp = parentSpGrp[0].upper() + parentSpGrp[1:].lower().rstrip('rh ')
2396#                 WriteCIFitem(self.fp, '_parent_space_group.name_H-M_alt',parentSpGrp)
2397# #                [Trans,Uvec,Vvec] = phasedict['General']['SGData']['fromParent']         #save these
2398#                 spacegroup = phasedict['General']['SGData']['MagSpGrp'].strip()
2399#                 spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
2400#                 WriteCIFitem(self.fp, '_space_group_magn.name_BNS',spacegroup)
2401#                 WriteCIFitem(self.fp, '_space_group.magn_point_group',phasedict['General']['SGData']['MagPtGp'])
2402
2403#                 # generate symmetry operations including centering and center of symmetry
2404#                 SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
2405#                     phasedict['General']['SGData'])
2406#                 SpnFlp = phasedict['General']['SGData']['SpnFlp']
2407#                 WriteCIFitem(self.fp, 'loop_\n    _space_group_symop_magn_operation.id\n    _space_group_symop_magn_operation.xyz')
2408#                 for i,op in enumerate(SymOpList,start=1):
2409#                     if SpnFlp[i-1] >0:
2410#                         opr = op.lower()+',+1'
2411#                     else:
2412#                         opr = op.lower()+',-1'
2413#                     WriteCIFitem(self.fp, '   {:3d}  {:}'.format(i,opr))
2414                   
2415            lam = None
2416            if 'X' in histblk['Instrument Parameters'][0]['Type'][0]:
2417                for k in ('Lam','Lam1'):
2418                    if k in histblk['Instrument Parameters'][0]:
2419                        lam = histblk['Instrument Parameters'][0][k][0]
2420                        break
2421            keV = None
2422            if lam: keV = 12.397639/lam
2423            # report cell contents                       
2424            WriteComposition(self.fp, self.Phases[phasenam], phasenam, self.parmDict, False, keV)
2425       
2426        def WriteSeqPhaseVals(phasenam,phasedict,pId,histname):
2427            'Write out the phase information for the selected phase'
2428            WriteCIFitem(self.fp, '_pd_phase_name', phasenam)
2429            cellList,cellSig = getCellwStrain(phasedict,self.seqData,pId,histname)
2430            T = self.Histograms[histname]['Sample Parameters']['Temperature']
2431            try:
2432                T = G2mth.ValEsd(T,-1.0)
2433            except:
2434                pass
2435            WriteCIFitem(self.fp, '_symmetry_cell_setting',
2436                         phasedict['General']['SGData']['SGSys'])
2437
2438            # generate symmetry operations including centering and center of symmetry
2439            # note that this would be better in WriteSeqOverallPhaseInfo() so there could
2440            # be only one copy per phase
2441            if phasedict['General']['Type'] in ['nuclear','macromolecular']:
2442                spacegroup = phasedict['General']['SGData']['SpGrp'].strip()
2443                # regularize capitalization and remove trailing H/R
2444                spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
2445                WriteCIFitem(self.fp, '_symmetry_space_group_name_H-M',spacegroup)
2446   
2447                # generate symmetry operations including centering and center of symmetry
2448                SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
2449                    phasedict['General']['SGData'])
2450                WriteCIFitem(self.fp, 'loop_\n    _space_group_symop_id\n    _space_group_symop_operation_xyz')
2451                for i,op in enumerate(SymOpList,start=1):
2452                    WriteCIFitem(self.fp, '   {:3d}  {:}'.format(i,op.lower()))
2453            elif phasedict['General']['Type'] == 'magnetic':
2454                parentSpGrp = phasedict['General']['SGData']['SpGrp'].strip()
2455                parentSpGrp = parentSpGrp[0].upper() + parentSpGrp[1:].lower().rstrip('rh ')
2456                WriteCIFitem(self.fp, '_parent_space_group.name_H-M_alt',parentSpGrp)
2457#                [Trans,Uvec,Vvec] = phasedict['General']['SGData']['fromParent']         #save these
2458                spacegroup = phasedict['General']['SGData']['MagSpGrp'].strip()
2459                spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
2460                WriteCIFitem(self.fp, '_space_group_magn.name_BNS',spacegroup)
2461                WriteCIFitem(self.fp, '_space_group.magn_point_group',phasedict['General']['SGData']['MagPtGp'])
2462
2463                # generate symmetry operations including centering and center of symmetry
2464                SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
2465                    phasedict['General']['SGData'])
2466                SpnFlp = phasedict['General']['SGData']['SpnFlp']
2467                WriteCIFitem(self.fp, 'loop_\n    _space_group_symop_magn_operation.id\n    _space_group_symop_magn_operation.xyz')
2468                for i,op in enumerate(SymOpList,start=1):
2469                    if SpnFlp[i-1] >0:
2470                        opr = op.lower()+',+1'
2471                    else:
2472                        opr = op.lower()+',-1'
2473                    WriteCIFitem(self.fp, '   {:3d}  {:}'.format(i,opr))
2474
2475            WriteCIFitem(self.fp,"_cell_measurement_temperature",T)
2476            defsigL = 3*[-0.00001] + 3*[-0.001] + [-0.01] # significance to use when no sigma
2477            prevsig = 0
2478            for lbl,defsig,val,sig in zip(cellNames,defsigL,cellList,cellSig):
2479                if sig:
2480                    txt = G2mth.ValEsd(val,sig)
2481                    prevsig = -sig # use this as the significance for next value
2482                else:
2483                    txt = G2mth.ValEsd(val,min(defsig,prevsig),True)
2484                WriteCIFitem(self.fp, '_cell_'+lbl,txt)
2485
2486            mass = G2mth.getMass(phasedict['General'])
2487            Volume = cellList[6]
2488            density = mass/(0.6022137*Volume)
2489            WriteCIFitem(self.fp, '_exptl_crystal_density_diffrn',
2490                    G2mth.ValEsd(density,-0.001))
2491
2492            # report atom params
2493            if phasedict['General']['Type'] in ['nuclear','macromolecular']:        #this needs macromolecular variant, etc!
2494                WriteSeqAtomsNuclear(self.fp, cellList, phasedict, phasenam, histname, 
2495                                         self.seqData, self.OverallParms['Rigid bodies'])
2496            else:
2497                print("Warning: no export for sequential "+str(phasedict['General']['Type'])+" coordinates implemented")
2498#                raise Exception("no export for "+str(phasedict['General']['Type'])+" coordinates implemented")
2499
2500            if phasedict['General']['Type'] == 'nuclear':
2501                WriteSeqDistances(phasenam,histname,phasedict,cellList,self.seqData)
2502
2503            # N.B. map info probably not possible w/sequential
2504#            if 'Map' in phasedict['General'] and 'minmax' in phasedict['General']['Map']:
2505#                WriteCIFitem(self.fp, '\n# Difference density results')
2506#                MinMax = phasedict['General']['Map']['minmax']
2507#                WriteCIFitem(self.fp, '_refine_diff_density_max',G2mth.ValEsd(MinMax[0],-0.009))
2508#                WriteCIFitem(self.fp, '_refine_diff_density_min',G2mth.ValEsd(MinMax[1],-0.009))
2509               
2510        def WritePhaseInfo(phasenam,quick=True,oneblock=True):
2511            'Write out the phase information for the selected phase'
2512            WriteCIFitem(self.fp, '\n# phase info for '+str(phasenam) + ' follows')
2513            phasedict = self.Phases[phasenam] # pointer to current phase info
2514            WriteCIFitem(self.fp, '_pd_phase_name', phasenam)
2515            cellList,cellSig = self.GetCell(phasenam)
2516            if quick:  # leave temperature as unknown
2517                WriteCIFitem(self.fp,"_cell_measurement_temperature","?")
2518            elif oneblock:
2519                pass # temperature should be written when the histogram saved later
2520            else: # get T set in _SelectPhaseT_CellSelectHist and possibly get new cell params
2521                T,hRanId = self.CellHistSelection.get(phasedict['ranId'],
2522                                                          ('?',None))
2523                try:
2524                    T = G2mth.ValEsd(T,-1.0)
2525                except:
2526                    pass
2527                WriteCIFitem(self.fp,"_cell_measurement_temperature",T)
2528                for h in self.Histograms:
2529                    if self.Histograms[h]['ranId'] == hRanId:
2530                        pId = phasedict['pId']
2531                        hId = self.Histograms[h]['hId']
2532                        cellList,cellSig = G2stIO.getCellSU(pId,hId,
2533                                        phasedict['General']['SGData'],
2534                                        self.parmDict,
2535                                        self.OverallParms['Covariance'])
2536                        break
2537                else:
2538                    T = '?'
2539
2540            defsigL = 3*[-0.00001] + 3*[-0.001] + [-0.01] # significance to use when no sigma
2541            prevsig = 0
2542            for lbl,defsig,val,sig in zip(cellNames,defsigL,cellList,cellSig):
2543                if sig:
2544                    txt = G2mth.ValEsd(val,sig)
2545                    prevsig = -sig # use this as the significance for next value
2546                else:
2547                    txt = G2mth.ValEsd(val,min(defsig,prevsig),True)
2548                WriteCIFitem(self.fp, '_cell_'+lbl,txt)
2549               
2550            density = G2mth.getDensity(phasedict['General'])[0]
2551            WriteCIFitem(self.fp, '_exptl_crystal_density_diffrn',
2552                    G2mth.ValEsd(density,-0.001))                   
2553
2554            WriteCIFitem(self.fp, '_symmetry_cell_setting',
2555                         phasedict['General']['SGData']['SGSys'])
2556
2557            if phasedict['General']['Type'] in ['nuclear','macromolecular']:
2558                spacegroup = phasedict['General']['SGData']['SpGrp'].strip()
2559                # regularize capitalization and remove trailing H/R
2560                spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
2561                WriteCIFitem(self.fp, '_symmetry_space_group_name_H-M',spacegroup)
2562   
2563                # generate symmetry operations including centering and center of symmetry
2564                SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
2565                    phasedict['General']['SGData'])
2566                WriteCIFitem(self.fp, 'loop_\n    _space_group_symop_id\n    _space_group_symop_operation_xyz')
2567                for i,op in enumerate(SymOpList,start=1):
2568                    WriteCIFitem(self.fp, '   {:3d}  {:}'.format(i,op.lower()))
2569            elif phasedict['General']['Type'] == 'magnetic':
2570                parentSpGrp = phasedict['General']['SGData']['SpGrp'].strip()
2571                parentSpGrp = parentSpGrp[0].upper() + parentSpGrp[1:].lower().rstrip('rh ')
2572                WriteCIFitem(self.fp, '_parent_space_group.name_H-M_alt',parentSpGrp)
2573#                [Trans,Uvec,Vvec] = phasedict['General']['SGData']['fromParent']         #save these
2574                spacegroup = phasedict['General']['SGData']['MagSpGrp'].strip()
2575                spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
2576                WriteCIFitem(self.fp, '_space_group_magn.name_BNS',spacegroup)
2577                WriteCIFitem(self.fp, '_space_group.magn_point_group',phasedict['General']['SGData']['MagPtGp'])
2578
2579                # generate symmetry operations including centering and center of symmetry
2580                SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
2581                    phasedict['General']['SGData'])
2582                SpnFlp = phasedict['General']['SGData']['SpnFlp']
2583                WriteCIFitem(self.fp, 'loop_\n    _space_group_symop_magn_operation.id\n    _space_group_symop_magn_operation.xyz')
2584                for i,op in enumerate(SymOpList,start=1):
2585                    if SpnFlp[i-1] >0:
2586                        opr = op.lower()+',+1'
2587                    else:
2588                        opr = op.lower()+',-1'
2589                    WriteCIFitem(self.fp, '   {:3d}  {:}'.format(i,opr))
2590
2591            # loop over histogram(s) used in this phase
2592            if not oneblock and not self.quickmode:
2593                # report pointers to the histograms used in this phase
2594                histlist = []
2595                for hist in self.Phases[phasenam]['Histograms']:
2596                    if self.Phases[phasenam]['Histograms'][hist]['Use']:
2597                        if phasebyhistDict.get(hist):
2598                            phasebyhistDict[hist].append(phasenam)
2599                        else:
2600                            phasebyhistDict[hist] = [phasenam,]
2601                        blockid = datablockidDict.get(hist)
2602                        if not blockid:
2603                            print("Internal error: no block for data. Phase "+str(
2604                                phasenam)+" histogram "+str(hist))
2605                            histlist = []
2606                            break
2607                        histlist.append(blockid)
2608
2609                if len(histlist) == 0:
2610                    WriteCIFitem(self.fp, '# Note: phase has no associated data')
2611
2612            # report atom params
2613            if phasedict['General']['Type'] in ['nuclear','macromolecular']:        #this needs macromolecular variant, etc!
2614                try:
2615                    self.labellist
2616                except AttributeError:
2617                    self.labellist = []
2618                WriteAtomsNuclear(self.fp, self.Phases[phasenam], phasenam,
2619                                  self.parmDict, self.sigDict, self.labellist,
2620                                      self.OverallParms['Rigid bodies'])
2621            else:
2622                try:
2623                    self.labellist
2624                except AttributeError:
2625                    self.labellist = []
2626                WriteAtomsMagnetic(self.fp, self.Phases[phasenam], phasenam,
2627                                  self.parmDict, self.sigDict, self.labellist)
2628#                raise Exception("no export for "+str(phasedict['General']['Type'])+" coordinates implemented")
2629            keV = None             
2630            if oneblock: # get xray wavelength
2631                lamlist = []
2632                for hist in self.Histograms:
2633                    if 'X' not in self.Histograms[hist]['Instrument Parameters'][0]['Type'][0]:
2634                        continue
2635                    for k in ('Lam','Lam1'):
2636                        if k in self.Histograms[hist]['Instrument Parameters'][0]:
2637                            lamlist.append(self.Histograms[hist]['Instrument Parameters'][0][k][0])
2638                            break
2639                if len(lamlist) == 1:
2640                    keV = 12.397639/lamlist[0]
2641
2642            # report cell contents                       
2643            WriteComposition(self.fp, self.Phases[phasenam], phasenam, self.parmDict, self.quickmode, keV)
2644            if not self.quickmode and phasedict['General']['Type'] == 'nuclear':      # report distances and angles
2645                WriteDistances(phasenam)
2646            if 'Map' in phasedict['General'] and 'minmax' in phasedict['General']['Map']:
2647                WriteCIFitem(self.fp, '\n# Difference density results')
2648                MinMax = phasedict['General']['Map']['minmax']
2649                WriteCIFitem(self.fp, '_refine_diff_density_max',G2mth.ValEsd(MinMax[0],-0.009))
2650                WriteCIFitem(self.fp, '_refine_diff_density_min',G2mth.ValEsd(MinMax[1],-0.009))
2651
2652        def WritePhaseInfoMM(phasenam,quick=True,oneblock=True):
2653            'Write out the phase information for the selected phase for a macromolecular phase'
2654            WriteCIFitem(self.fp, '\n# phase info for '+str(phasenam) + ' follows')
2655            phasedict = self.Phases[phasenam] # pointer to current phase info
2656            WriteCIFitem(self.fp, '_cell.entry_id', phasenam)
2657            cellList,cellSig = self.GetCell(phasenam)
2658            if oneblock:
2659                pass # temperature should be written when the histogram saved later
2660            else: # get T set in _SelectPhaseT_CellSelectHist and possibly get new cell params
2661                T,hRanId = self.CellHistSelection.get(phasedict['ranId'],
2662                                                          ('?',None))
2663                try:
2664                    T = G2mth.ValEsd(T,-1.0)
2665                except:
2666                    pass
2667                WriteCIFitem(self.fp,"_cell_measurement.temp",T)
2668                for h in self.Histograms:
2669                    if self.Histograms[h]['ranId'] == hRanId:
2670                        pId = phasedict['pId']
2671                        hId = self.Histograms[h]['hId']
2672                        cellList,cellSig = G2stIO.getCellSU(pId,hId,
2673                                        phasedict['General']['SGData'],
2674                                        self.parmDict,
2675                                        self.OverallParms['Covariance'])
2676                        break
2677                else:
2678                    T = '?'
2679
2680            defsigL = 3*[-0.00001] + 3*[-0.001] + [-0.01] # significance to use when no sigma
2681            prevsig = 0
2682            for lbl,defsig,val,sig in zip(cellNames,defsigL,cellList,cellSig):
2683                if sig:
2684                    txt = G2mth.ValEsd(val,sig)
2685                    prevsig = -sig # use this as the significance for next value
2686                else:
2687                    txt = G2mth.ValEsd(val,min(defsig,prevsig),True)
2688                WriteCIFitem(self.fp, '_cell.'+lbl,txt)
2689               
2690            density = G2mth.getDensity(phasedict['General'])[0]
2691            WriteCIFitem(self.fp, '_exptl_crystal.density_diffrn',
2692                    G2mth.ValEsd(density,-0.001))                   
2693
2694            WriteCIFitem(self.fp, '_symmetry.cell_setting',
2695                         phasedict['General']['SGData']['SGSys'])
2696
2697            spacegroup = phasedict['General']['SGData']['SpGrp'].strip()
2698            # regularize capitalization and remove trailing H/R
2699            spacegroup = spacegroup[0].upper() + spacegroup[1:].lower().rstrip('rh ')
2700            WriteCIFitem(self.fp, '_symmetry.space_group_name_H-M',spacegroup)
2701
2702            # generate symmetry operations including centering and center of symmetry
2703            SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
2704                phasedict['General']['SGData'])
2705            WriteCIFitem(self.fp, 'loop_\n    _space_group.symop_id\n    _space_group.symop_operation_xyz')
2706            for i,op in enumerate(SymOpList,start=1):
2707                WriteCIFitem(self.fp, '   {:3d}  {:}'.format(i,op.lower()))
2708
2709            # loop over histogram(s) used in this phase
2710            if not oneblock and not self.quickmode:
2711                # report pointers to the histograms used in this phase
2712                histlist = []
2713                for hist in self.Phases[phasenam]['Histograms']:
2714                    if self.Phases[phasenam]['Histograms'][hist]['Use']:
2715                        if phasebyhistDict.get(hist):
2716                            phasebyhistDict[hist].append(phasenam)
2717                        else:
2718                            phasebyhistDict[hist] = [phasenam,]
2719                        blockid = datablockidDict.get(hist)
2720                        if not blockid:
2721                            print("Internal error: no block for data. Phase "+str(
2722                                phasenam)+" histogram "+str(hist))
2723                            histlist = []
2724                            break
2725                        histlist.append(blockid)
2726
2727                if len(histlist) == 0:
2728                    WriteCIFitem(self.fp, '# Note: phase has no associated data')
2729
2730            # report atom params
2731            try:
2732                self.labellist
2733            except AttributeError:
2734                self.labellist = []
2735            WriteAtomsMM(self.fp, self.Phases[phasenam], phasenam,
2736                              self.parmDict, self.sigDict, 
2737                                  self.OverallParms['Rigid bodies'])
2738
2739            keV = None             
2740            if oneblock: # get xray wavelength
2741                lamlist = []
2742                for hist in self.Histograms:
2743                    if 'X' not in self.Histograms[hist]['Instrument Parameters'][0]['Type'][0]:
2744                        continue
2745                    for k in ('Lam','Lam1'):
2746                        if k in self.Histograms[hist]['Instrument Parameters'][0]:
2747                            lamlist.append(self.Histograms[hist]['Instrument Parameters'][0][k][0])
2748                            break
2749                if len(lamlist) == 1:
2750                    keV = 12.397639/lamlist[0]
2751
2752            # report cell contents                       
2753            WriteCompositionMM(self.fp, self.Phases[phasenam], phasenam, self.parmDict, self.quickmode, keV)
2754            #if not self.quickmode and phasedict['General']['Type'] == 'nuclear':      # report distances and angles
2755            #    WriteDistances(phasenam)
2756            if 'Map' in phasedict['General'] and 'minmax' in phasedict['General']['Map']:
2757                WriteCIFitem(self.fp, '\n# Difference density results')
2758                MinMax = phasedict['General']['Map']['minmax']
2759                WriteCIFitem(self.fp, '_refine.diff_density_max',G2mth.ValEsd(MinMax[0],-0.009))
2760                WriteCIFitem(self.fp, '_refine.diff_density_min',G2mth.ValEsd(MinMax[1],-0.009))
2761               
2762        def Yfmt(ndec,val):
2763            'Format intensity values'
2764            try:
2765                out = ("{:."+str(ndec)+"f}").format(val)
2766                out = out.rstrip('0')  # strip zeros to right of decimal
2767                return out.rstrip('.')  # and decimal place when not needed
2768            except TypeError:
2769                print(val)
2770                return '.'
2771           
2772        def WriteReflStat(refcount,hklmin,hklmax,dmin,dmax,nRefSets=1):
2773            'Write reflection statistics'
2774            WriteCIFitem(self.fp, '_reflns_number_total', str(refcount))
2775            if hklmin is not None and nRefSets == 1: # hkl range has no meaning with multiple phases
2776                WriteCIFitem(self.fp, '_reflns_limit_h_min', str(int(hklmin[0])))
2777                WriteCIFitem(self.fp, '_reflns_limit_h_max', str(int(hklmax[0])))
2778                WriteCIFitem(self.fp, '_reflns_limit_k_min', str(int(hklmin[1])))
2779                WriteCIFitem(self.fp, '_reflns_limit_k_max', str(int(hklmax[1])))
2780                WriteCIFitem(self.fp, '_reflns_limit_l_min', str(int(hklmin[2])))
2781                WriteCIFitem(self.fp, '_reflns_limit_l_max', str(int(hklmax[2])))
2782            if hklmin is not None:
2783                WriteCIFitem(self.fp, '_reflns_d_resolution_low  ', G2mth.ValEsd(dmax,-0.009))
2784                WriteCIFitem(self.fp, '_reflns_d_resolution_high ', G2mth.ValEsd(dmin,-0.009))
2785
2786        def WritePowderData(histlbl,seq=False):
2787            'Write out the selected powder diffraction histogram info'
2788            histblk = self.Histograms[histlbl]
2789            inst = histblk['Instrument Parameters'][0]
2790            if seq:
2791                resdblk = histblk['Residuals']
2792            else:
2793                resdblk = histblk
2794            hId = histblk['hId']
2795            pfx = ':' + str(hId) + ':'
2796
2797            if 'Lam1' in inst:
2798                ratio = self.parmDict.get('I(L2)/I(L1)',inst['I(L2)/I(L1)'][1])
2799                sratio = self.sigDict.get('I(L2)/I(L1)',-0.0009)
2800                lam1 = self.parmDict.get('Lam1',inst['Lam1'][1])
2801                slam1 = self.sigDict.get('Lam1',-0.00009)
2802                lam2 = self.parmDict.get('Lam2',inst['Lam2'][1])
2803                slam2 = self.sigDict.get('Lam2',-0.00009)
2804                # always assume Ka1 & Ka2 if two wavelengths are present
2805                WriteCIFitem(self.fp, '_diffrn_radiation_type','K\\a~1,2~')
2806                WriteCIFitem(self.fp, 'loop_' +
2807                             '\n   _diffrn_radiation_wavelength' +
2808                             '\n   _diffrn_radiation_wavelength_wt' +
2809                             '\n   _diffrn_radiation_wavelength_id')
2810                WriteCIFitem(self.fp, '  ' + PutInCol(G2mth.ValEsd(lam1,slam1),15)+
2811                             PutInCol('1.0',15) +
2812                             PutInCol('1',5))
2813                WriteCIFitem(self.fp, '  ' + PutInCol(G2mth.ValEsd(lam2,slam2),15)+
2814                             PutInCol(G2mth.ValEsd(ratio,sratio),15)+
2815                             PutInCol('2',5))
2816            elif 'Lam' in inst:
2817                lam1 = self.parmDict.get('Lam',inst['Lam'][1])
2818                slam1 = self.sigDict.get('Lam',-0.00009)
2819                WriteCIFitem(self.fp, '_diffrn_radiation_wavelength',G2mth.ValEsd(lam1,slam1))
2820
2821            if not oneblock:
2822                if not phasebyhistDict.get(histlbl) and not seq:
2823                    WriteCIFitem(self.fp, '\n# No phases associated with this data set')
2824                elif len(self.Phases) == 1:
2825                    pId = self.Phases[list(self.Phases.keys())[0]]['pId']
2826                    pfx = str(pId)+':'+str(hId)+':'
2827                    WriteCIFitem(self.fp, '_refine_ls_R_F_factor      ','%.5f'%(resdblk[pfx+'Rf']/100.))
2828                    WriteCIFitem(self.fp, '_refine_ls_R_Fsqd_factor   ','%.5f'%(resdblk[pfx+'Rf^2']/100.))
2829                else:
2830                    WriteCIFitem(self.fp, '\n# PHASE TABLE')
2831                    WriteCIFitem(self.fp, 'loop_' +
2832                                 '\n   _pd_phase_id' +
2833                                 '\n   _pd_phase_block_id' +
2834                                 '\n   _pd_phase_mass_%')
2835                    wtFrSum = 0.
2836                    for phasenam in phasebyhistDict.get(histlbl):
2837                        hapData = self.Phases[phasenam]['Histograms'][histlbl]
2838                        General = self.Phases[phasenam]['General']
2839                        wtFrSum += hapData['Scale'][0]*General['Mass']
2840
2841                    for phasenam in phasebyhistDict.get(histlbl):
2842                        hapData = self.Phases[phasenam]['Histograms'][histlbl]
2843                        General = self.Phases[phasenam]['General']
2844                        wtFr = hapData['Scale'][0]*General['Mass']/wtFrSum
2845                        pfx = str(self.Phases[phasenam]['pId'])+':'+str(hId)+':'
2846                        if pfx+'Scale' in self.sigDict:
2847                            sig = self.sigDict[pfx+'Scale']*wtFr/hapData['Scale'][0]
2848                        else:
2849                            sig = -0.0001
2850                        WriteCIFitem(self.fp,
2851                            '  '+
2852                            str(self.Phases[phasenam]['pId']) +
2853                            '  '+datablockidDict[phasenam]+
2854                            '  '+G2mth.ValEsd(wtFr,sig)
2855                            )
2856                    WriteCIFitem(self.fp, 'loop_' +
2857                                 '\n   _gsas_proc_phase_R_F_factor' +
2858                                 '\n   _gsas_proc_phase_R_Fsqd_factor' +
2859                                 '\n   _gsas_proc_phase_id' +
2860                                 '\n   _gsas_proc_phase_block_id')
2861                    for phasenam in phasebyhistDict.get(histlbl):
2862                        pfx = str(self.Phases[phasenam]['pId'])+':'+str(hId)+':'
2863                        WriteCIFitem(self.fp,
2864                            '  '+
2865                            '  '+G2mth.ValEsd(resdblk[pfx+'Rf']/100.,-.00009) +
2866                            '  '+G2mth.ValEsd(resdblk[pfx+'Rf^2']/100.,-.00009)+
2867                            '  '+str(self.Phases[phasenam]['pId'])+
2868                            '  '+datablockidDict[phasenam]
2869                            )
2870            elif len(self.Phases) == 1:
2871                # single phase in this histogram
2872                # get the phase number here
2873                pId = self.Phases[list(self.Phases.keys())[0]]['pId']
2874                pfx = str(pId)+':'+str(hId)+':'
2875                WriteCIFitem(self.fp, '_refine_ls_R_F_factor      ','%.5f'%(resdblk[pfx+'Rf']/100.))
2876                WriteCIFitem(self.fp, '_refine_ls_R_Fsqd_factor   ','%.5f'%(resdblk[pfx+'Rf^2']/100.))
2877               
2878            try:
2879                WriteCIFitem(self.fp, '_pd_proc_ls_prof_R_factor   ','%.5f'%(resdblk['R']/100.))
2880                WriteCIFitem(self.fp, '_pd_proc_ls_prof_wR_factor  ','%.5f'%(resdblk['wR']/100.))
2881                WriteCIFitem(self.fp, '_gsas_proc_ls_prof_R_B_factor ','%.5f'%(resdblk['Rb']/100.))
2882                WriteCIFitem(self.fp, '_gsas_proc_ls_prof_wR_B_factor','%.5f'%(resdblk['wRb']/100.))
2883                WriteCIFitem(self.fp, '_pd_proc_ls_prof_wR_expected','%.5f'%(resdblk['wRmin']/100.))
2884                if not oneblock: # written in WriteOverall, don't repeat in a one-block CIF
2885                    WriteCIFitem(self.fp, '_refine_ls_goodness_of_fit_all','%.2f'%(resdblk['wR']/resdblk['wRmin']))
2886            except KeyError:
2887                pass
2888
2889            if histblk['Instrument Parameters'][0]['Type'][1][1] == 'X':
2890                WriteCIFitem(self.fp, '_diffrn_radiation_probe','x-ray')
2891                pola = histblk['Instrument Parameters'][0].get('Polariz.')
2892                if pola:
2893                    pfx = ':' + str(hId) + ':'
2894                    sig = self.sigDict.get(pfx+'Polariz.',-0.0009)
2895                    txt = G2mth.ValEsd(pola[1],sig)
2896                    WriteCIFitem(self.fp, '_diffrn_radiation_polarisn_ratio',txt)
2897            elif histblk['Instrument Parameters'][0]['Type'][1][1] == 'N':
2898                WriteCIFitem(self.fp, '_diffrn_radiation_probe','neutron')
2899            if 'T' in inst['Type'][0]:
2900                txt = G2mth.ValEsd(inst['2-theta'][0],-0.009)
2901                WriteCIFitem(self.fp, '_pd_meas_2theta_fixed',txt)
2902
2903            WriteCIFitem(self.fp, '_pd_proc_ls_background_function',FormatBackground(histblk['Background'],histblk['hId']))
2904
2905            # TODO: this will need help from Bob
2906            #WriteCIFitem(self.fp, '_exptl_absorpt_process_details','?')
2907            #WriteCIFitem(self.fp, '_exptl_absorpt_correction_T_min','?')
2908            #WriteCIFitem(self.fp, '_exptl_absorpt_correction_T_max','?')
2909            #C extinction
2910            #WRITE(IUCIF,'(A)') '# Extinction correction'
2911            #CALL WRVAL(IUCIF,'_gsas_exptl_extinct_corr_T_min',TEXT(1:10))
2912            #CALL WRVAL(IUCIF,'_gsas_exptl_extinct_corr_T_max',TEXT(11:20))
2913
2914            # code removed because it is causing duplication in histogram block 1/26/19 BHT
2915            #if not oneblock:                 # instrumental profile terms go here
2916            #    WriteCIFitem(self.fp, '_pd_proc_ls_profile_function',
2917            #        FormatInstProfile(histblk["Instrument Parameters"],histblk['hId']))
2918
2919            #refprx = '_refln.' # mm
2920            refprx = '_refln_' # normal
2921            # data collection parameters for the powder dataset
2922
2923            temperature = histblk['Sample Parameters'].get('Temperature') # G2 uses K
2924            if not temperature:
2925                T = '?'
2926            else:
2927                T = G2mth.ValEsd(temperature,-0.009,True) # CIF uses K
2928            WriteCIFitem(self.fp, '_diffrn_ambient_temperature',T)
2929
2930            pressure = histblk['Sample Parameters'].get('Pressure') #G2 uses mega-Pascal
2931            if not pressure:
2932                P = '?'
2933            else:
2934                P = G2mth.ValEsd(pressure*1000,-0.09,True) # CIF uses kilopascal (G2 Mpa)
2935            WriteCIFitem(self.fp, '_diffrn_ambient_pressure',P)
2936
2937            WriteCIFitem(self.fp, '\n# STRUCTURE FACTOR TABLE')
2938            # compute maximum intensity reflection
2939            Imax = 0
2940            phaselist = []
2941            for phasenam in histblk['Reflection Lists']:
2942                try:
2943                    scale = self.Phases[phasenam]['Histograms'][histlbl]['Scale'][0]
2944                except KeyError: # reflection table from removed phase?
2945                    continue
2946                phaselist.append(phasenam)
2947                refList = np.asarray(histblk['Reflection Lists'][phasenam]['RefList'])
2948                I100 = scale*refList.T[8]*refList.T[11]
2949                #Icorr = np.array([refl[13] for refl in histblk['Reflection Lists'][phasenam]])[0]
2950                #FO2 = np.array([refl[8] for refl in histblk['Reflection Lists'][phasenam]])
2951                #I100 = scale*FO2*Icorr
2952                Imax = max(Imax,max(I100))
2953
2954            WriteCIFitem(self.fp, 'loop_')
2955            if len(phaselist) > 1:
2956                WriteCIFitem(self.fp, '   _pd_refln_phase_id')
2957            WriteCIFitem(self.fp, '   ' + refprx + 'index_h' +
2958                         '\n   ' + refprx + 'index_k' +
2959                         '\n   ' + refprx + 'index_l' +
2960                         '\n   ' + refprx + 'F_squared_meas' +
2961                         '\n   ' + refprx + 'F_squared_calc' +
2962                         '\n   ' + refprx + 'phase_calc' +
2963                         '\n   _refln_d_spacing')
2964            if Imax > 0:
2965                WriteCIFitem(self.fp, '   _gsas_i100_meas')
2966
2967            refcount = 0
2968            hklmin = None
2969            hklmax = None
2970            dmax = None
2971            dmin = None
2972            for phasenam in phaselist:
2973                scale = self.Phases[phasenam]['Histograms'][histlbl]['Scale'][0]
2974                phaseid = self.Phases[phasenam]['pId']
2975                refcount += len(histblk['Reflection Lists'][phasenam]['RefList'])
2976                refList = np.asarray(histblk['Reflection Lists'][phasenam]['RefList'])
2977                I100 = scale*refList.T[8]*refList.T[11]
2978                for j,ref in enumerate(histblk['Reflection Lists'][phasenam]['RefList']):
2979                    if DEBUG:
2980                        print('DEBUG: skipping reflection list')
2981                        break
2982                    if hklmin is None:
2983                        hklmin = copy.copy(ref[0:3])
2984                        hklmax = copy.copy(ref[0:3])
2985                    if dmin is None:
2986                         dmax = dmin = ref[4]
2987                    if len(phaselist) > 1:
2988                        s = PutInCol(phaseid,2)
2989                    else:
2990                        s = ""
2991                    for i,hkl in enumerate(ref[0:3]):
2992                        hklmax[i] = max(hkl,hklmax[i])
2993                        hklmin[i] = min(hkl,hklmin[i])
2994                        s += PutInCol(int(hkl),4)
2995                    for I in ref[8:10]:
2996                        s += PutInCol(G2mth.ValEsd(I,-0.0009),10)
2997                    s += PutInCol(G2mth.ValEsd(ref[10],-0.9),7)
2998                    dmax = max(dmax,ref[4])
2999                    dmin = min(dmin,ref[4])
3000                    s += PutInCol(G2mth.ValEsd(ref[4],-0.00009),8)
3001                    if Imax > 0:
3002                        s += PutInCol(G2mth.ValEsd(100.*I100[j]/Imax,-0.09),6)
3003                    WriteCIFitem(self.fp, "  "+s)
3004
3005            WriteReflStat(refcount,hklmin,hklmax,dmin,dmax,len(phaselist))
3006            WriteCIFitem(self.fp, '\n# POWDER DATA TABLE')
3007            # is data fixed step? If the step varies by <0.01% treat as fixed step
3008            steps = abs(histblk['Data'][0][1:] - histblk['Data'][0][:-1])
3009            if (max(steps)-min(steps)) > np.mean(steps)/10000.:
3010                fixedstep = False
3011            else:
3012                fixedstep = True
3013
3014            zero = None
3015            if fixedstep and 'T' not in inst['Type'][0]: # and not TOF
3016                WriteCIFitem(self.fp, '_pd_meas_2theta_range_min', G2mth.ValEsd(histblk['Data'][0][0],-0.00009))
3017                WriteCIFitem(self.fp, '_pd_meas_2theta_range_max', G2mth.ValEsd(histblk['Data'][0][-1],-0.00009))
3018                WriteCIFitem(self.fp, '_pd_meas_2theta_range_inc', G2mth.ValEsd(np.mean(steps),-0.00009))
3019                # zero correct, if defined
3020                zerolst = histblk['Instrument Parameters'][0].get('Zero')
3021                if zerolst: zero = zerolst[1]
3022                zero = self.parmDict.get('Zero',zero)
3023                if zero:
3024                    WriteCIFitem(self.fp, '_pd_proc_2theta_range_min', G2mth.ValEsd(histblk['Data'][0][0]-zero,-0.00009))
3025                    WriteCIFitem(self.fp, '_pd_proc_2theta_range_max', G2mth.ValEsd(histblk['Data'][0][-1]-zero,-0.00009))
3026                    WriteCIFitem(self.fp, '_pd_proc_2theta_range_inc', G2mth.ValEsd(steps.sum()/len(steps),-0.00009))
3027
3028            if zero:
3029                WriteCIFitem(self.fp, '_pd_proc_number_of_points', str(len(histblk['Data'][0])))
3030            else:
3031                WriteCIFitem(self.fp, '_pd_meas_number_of_points', str(len(histblk['Data'][0])))
3032            WriteCIFitem(self.fp, '\nloop_')
3033            #            WriteCIFitem(self.fp, '   _pd_proc_d_spacing') # need easy way to get this
3034            if not fixedstep:
3035                if zero:
3036                    WriteCIFitem(self.fp, '   _pd_proc_2theta_corrected')
3037                elif 'T' in inst['Type'][0]: # and not TOF
3038                    WriteCIFitem(self.fp, '   _pd_meas_time_of_flight')
3039                else:
3040                    WriteCIFitem(self.fp, '   _pd_meas_2theta_scan')
3041            # at least for now, always report weights.
3042            #if countsdata:
3043            #    WriteCIFitem(self.fp, '   _pd_meas_counts_total')
3044            #else:
3045            WriteCIFitem(self.fp, '   _pd_meas_intensity_total')
3046            WriteCIFitem(self.fp, '   _pd_calc_intensity_total')
3047            WriteCIFitem(self.fp, '   _pd_proc_intensity_bkg_calc')
3048            WriteCIFitem(self.fp, '   _pd_proc_ls_weight')
3049            maxY = max(histblk['Data'][1].max(),histblk['Data'][3].max())
3050            if maxY < 0: maxY *= -10 # this should never happen, but...
3051            ndec = max(0,10-int(np.log10(maxY))-1) # 10 sig figs should be enough
3052            maxSU = histblk['Data'][2].max()
3053            if maxSU < 0: maxSU *= -1 # this should never happen, but...
3054            ndecSU = max(0,8-int(np.log10(maxSU))-1) # 8 sig figs should be enough
3055            lowlim,highlim = histblk['Limits'][1]
3056
3057            if DEBUG:
3058                print('DEBUG: skipping profile list')
3059            else:
3060                for x,yobs,yw,ycalc,ybkg in zip(histblk['Data'][0].data,        #get the data from these masked arrays
3061                                                histblk['Data'][1].data,
3062                                                histblk['Data'][2].data,
3063                                                histblk['Data'][3].data,
3064                                                histblk['Data'][4].data):
3065                    if lowlim <= x <= highlim:
3066                        pass
3067                    else:
3068                        yw = 0.0 # show the point is not in use
3069
3070                    if fixedstep:
3071                        s = ""
3072                    elif zero:
3073                        s = PutInCol(G2mth.ValEsd(x-zero,-0.00009),10)
3074                    else:
3075                        s = PutInCol(G2mth.ValEsd(x,-0.00009),10)
3076                    s += PutInCol(Yfmt(ndec,yobs),12)
3077                    s += PutInCol(Yfmt(ndec,ycalc),12)
3078                    s += PutInCol(Yfmt(ndec,ybkg),11)
3079                    s += PutInCol(Yfmt(ndecSU,yw),9)
3080                    WriteCIFitem(self.fp, "  "+s)
3081                   
3082        def WritePowderDataMM(histlbl,seq=False):
3083            'Write out the selected powder diffraction histogram info'
3084            histblk = self.Histograms[histlbl]
3085            inst = histblk['Instrument Parameters'][0]
3086            hId = histblk['hId']
3087            pfx = ':' + str(hId) + ':'
3088
3089            WriteCIFitem(self.fp, '_diffrn.id',str(hId))
3090            WriteCIFitem(self.fp, '_diffrn.crystal_id',str(hId))
3091           
3092            if 'Lam1' in inst:
3093                ratio = self.parmDict.get('I(L2)/I(L1)',inst['I(L2)/I(L1)'][1])
3094                sratio = self.sigDict.get('I(L2)/I(L1)',-0.0009)
3095                lam1 = self.parmDict.get('Lam1',inst['Lam1'][1])
3096                slam1 = self.sigDict.get('Lam1',-0.00009)
3097                lam2 = self.parmDict.get('Lam2',inst['Lam2'][1])
3098                slam2 = self.sigDict.get('Lam2',-0.00009)
3099                # always assume Ka1 & Ka2 if two wavelengths are present
3100                WriteCIFitem(self.fp, '_diffrn_radiation.type','K\\a~1,2~')
3101                WriteCIFitem(self.fp, 'loop_' +
3102                             '\n   _diffrn_radiation_wavelength.wavelength' +
3103                             '\n   _diffrn_radiation_wavelength.wt' +
3104                             '\n   _diffrn_radiation_wavelength.id')
3105                WriteCIFitem(self.fp, '  ' + PutInCol(G2mth.ValEsd(lam1,slam1),15)+
3106                             PutInCol('1.0',15) +
3107                             PutInCol('1',5))
3108                WriteCIFitem(self.fp, '  ' + PutInCol(G2mth.ValEsd(lam2,slam2),15)+
3109                             PutInCol(G2mth.ValEsd(ratio,sratio),15)+
3110                             PutInCol('2',5))
3111            elif 'Lam' in inst:
3112                WriteCIFitem(self.fp, '_diffrn_radiation.diffrn_id',str(hId))
3113                WriteCIFitem(self.fp, '_diffrn_radiation.wavelength_id','1')
3114                WriteCIFitem(self.fp, '_diffrn_radiation_wavelength.id','1') 
3115                lam1 = self.parmDict.get('Lam',inst['Lam'][1])
3116                slam1 = self.sigDict.get('Lam',-0.00009)
3117                WriteCIFitem(self.fp, '_diffrn_radiation_wavelength.wavelength',G2mth.ValEsd(lam1,slam1))
3118
3119            if not oneblock:
3120                if seq:
3121                    pass
3122                elif not phasebyhistDict.get(histlbl):
3123                    WriteCIFitem(self.fp, '\n# No phases associated with this data set')
3124                else:
3125                    WriteCIFitem(self.fp, '\n# PHASE TABLE')
3126                    WriteCIFitem(self.fp, 'loop_' +
3127                                 '\n   _pd_phase_id' +
3128                                 '\n   _pd_phase_block_id' +
3129                                 '\n   _pd_phase_mass_%')
3130                    wtFrSum = 0.
3131                    for phasenam in phasebyhistDict.get(histlbl):
3132                        hapData = self.Phases[phasenam]['Histograms'][histlbl]
3133                        General = self.Phases[phasenam]['General']
3134                        wtFrSum += hapData['Scale'][0]*General['Mass']
3135
3136                    for phasenam in phasebyhistDict.get(histlbl):
3137                        hapData = self.Phases[phasenam]['Histograms'][histlbl]
3138                        General = self.Phases[phasenam]['General']
3139                        wtFr = hapData['Scale'][0]*General['Mass']/wtFrSum
3140                        pfx = str(self.Phases[phasenam]['pId'])+':'+str(hId)+':'
3141                        if pfx+'Scale' in self.sigDict:
3142                            sig = self.sigDict[pfx+'Scale']*wtFr/hapData['Scale'][0]
3143                        else:
3144                            sig = -0.0001
3145                        WriteCIFitem(self.fp,
3146                            '  '+
3147                            str(self.Phases[phasenam]['pId']) +
3148                            '  '+datablockidDict[phasenam]+
3149                            '  '+G2mth.ValEsd(wtFr,sig)
3150                            )
3151                    WriteCIFitem(self.fp, 'loop_' +
3152                                 '\n   _gsas_proc_phase_R_F_factor' +
3153                                 '\n   _gsas_proc_phase_R_Fsqd_factor' +
3154                                 '\n   _gsas_proc_phase_id' +
3155                                 '\n   _gsas_proc_phase_block_id')
3156                    for phasenam in phasebyhistDict.get(histlbl):
3157                        pfx = str(self.Phases[phasenam]['pId'])+':'+str(hId)+':'
3158                        WriteCIFitem(self.fp,
3159                            '  '+
3160                            '  '+G2mth.ValEsd(histblk[pfx+'Rf']/100.,-.00009) +
3161                            '  '+G2mth.ValEsd(histblk[pfx+'Rf^2']/100.,-.00009)+
3162                            '  '+str(self.Phases[phasenam]['pId'])+
3163                            '  '+datablockidDict[phasenam]
3164                            )
3165            elif len(self.Phases) == 1:
3166                # single phase in this histogram
3167                # get the phase number here
3168                pId = self.Phases[list(self.Phases.keys())[0]]['pId']
3169                pfx = str(pId)+':'+str(hId)+':'
3170                WriteCIFitem(self.fp, '_refine.ls_R_factor_all    ','%.5f'%(histblk[pfx+'Rf']/100.))
3171                WriteCIFitem(self.fp, '_refine_ls.R_Fsqd_factor   ','%.5f'%(histblk[pfx+'Rf^2']/100.))
3172               
3173            try:
3174                WriteCIFitem(self.fp, '_pd_proc_ls_prof_R_factor   ','%.5f'%(histblk['R']/100.))
3175                WriteCIFitem(self.fp, '_pd_proc_ls_prof_wR_factor  ','%.5f'%(histblk['wR']/100.))
3176                WriteCIFitem(self.fp, '_gsas_proc_ls_prof_R_B_factor ','%.5f'%(histblk['Rb']/100.))
3177                WriteCIFitem(self.fp, '_gsas_proc_ls_prof_wR_B_factor','%.5f'%(histblk['wRb']/100.))
3178                WriteCIFitem(self.fp, '_pd_proc_ls_prof_wR_expected','%.5f'%(histblk['wRmin']/100.))
3179            except KeyError:
3180                pass
3181
3182            if histblk['Instrument Parameters'][0]['Type'][1][1] == 'X':
3183                WriteCIFitem(self.fp, '_diffrn_radiation.probe','x-ray')
3184                pola = histblk['Instrument Parameters'][0].get('Polariz.')
3185                if pola:
3186                    pfx = ':' + str(hId) + ':'
3187                    sig = self.sigDict.get(pfx+'Polariz.',-0.0009)
3188                    txt = G2mth.ValEsd(pola[1],sig)
3189                    WriteCIFitem(self.fp, '_diffrn_radiation.polarisn_ratio',txt)
3190            elif histblk['Instrument Parameters'][0]['Type'][1][1] == 'N':
3191                WriteCIFitem(self.fp, '_diffrn_radiation.probe','neutron')
3192            if 'T' in inst['Type'][0]:
3193                txt = G2mth.ValEsd(inst['2-theta'][0],-0.009)
3194                WriteCIFitem(self.fp, '_pd_meas_2theta_fixed',txt)
3195
3196            WriteCIFitem(self.fp, '_pd_proc_ls_background_function',FormatBackground(histblk['Background'],histblk['hId']))
3197
3198            # TODO: this will need help from Bob
3199            #WriteCIFitem(self.fp, '_exptl_absorpt_process_details','?')
3200            #WriteCIFitem(self.fp, '_exptl_absorpt_correction_T_min','?')
3201            #WriteCIFitem(self.fp, '_exptl_absorpt_correction_T_max','?')
3202            #C extinction
3203            #WRITE(IUCIF,'(A)') '# Extinction correction'
3204            #CALL WRVAL(IUCIF,'_gsas_exptl_extinct_corr_T_min',TEXT(1:10))
3205            #CALL WRVAL(IUCIF,'_gsas_exptl_extinct_corr_T_max',TEXT(11:20))
3206
3207            # code removed because it is causing duplication in histogram block 1/26/19 BHT
3208            #if not oneblock:                 # instrumental profile terms go here
3209            #    WriteCIFitem(self.fp, '_pd_proc_ls_profile_function',
3210            #        FormatInstProfile(histblk["Instrument Parameters"],histblk['hId']))
3211
3212            # data collection parameters for the powder dataset
3213
3214            temperature = histblk['Sample Parameters'].get('Temperature') # G2 uses K
3215            if not temperature:
3216                T = '?'
3217            else:
3218                T = G2mth.ValEsd(temperature,-0.009,True) # CIF uses K
3219            WriteCIFitem(self.fp, '_diffrn_ambient.temp',T)
3220
3221            pressure = histblk['Sample Parameters'].get('Pressure') #G2 uses mega-Pascal
3222            if not pressure:
3223                P = '?'
3224            else:
3225                P = G2mth.ValEsd(pressure*1000,-0.09,True) # CIF uses kilopascal (G2 Mpa)
3226            WriteCIFitem(self.fp, '_diffrn_ambient.pressure',P)
3227
3228            WriteCIFitem(self.fp, '\n# STRUCTURE FACTOR TABLE')
3229            # compute maximum intensity reflection
3230            Imax = 0
3231            phaselist = []
3232            for phasenam in histblk['Reflection Lists']:
3233                try:
3234                    scale = self.Phases[phasenam]['Histograms'][histlbl]['Scale'][0]
3235                except KeyError: # reflection table from removed phase?
3236                    continue
3237                phaselist.append(phasenam)
3238                refList = np.asarray(histblk['Reflection Lists'][phasenam]['RefList'])
3239                I100 = scale*refList.T[8]*refList.T[11]
3240                #Icorr = np.array([refl[13] for refl in histblk['Reflection Lists'][phasenam]])[0]
3241                #FO2 = np.array([refl[8] for refl in histblk['Reflection Lists'][phasenam]])
3242                #I100 = scale*FO2*Icorr
3243#                Imax = max(Imax,max(I100))
3244
3245            WriteCIFitem(self.fp, 'loop_')
3246            #refprx = '_refln.' # mm
3247            if len(phaselist) > 1:
3248                WriteCIFitem(self.fp, '   _pd_refln_phase_id')
3249            WriteCIFitem(self.fp, '   _refln.index_h' +
3250                         '\n   _refln.index_k' +
3251                         '\n   _refln.index_l' +
3252                         '\n   _refln.F_squared_meas' +
3253                         '\n   _refln.F_squared_calc' +
3254                         '\n   _refln.phase_calc' +
3255                         '\n   _refln.d_spacing' +
3256                         '\n   _refln.status' +
3257                         '\n   _refln.crystal_id' +
3258                         '\n   _refln.wavelength_id' + 
3259                         '\n   _refln.scale_group_code' + 
3260                         '\n   _refln.F_squared_sigma')
3261           
3262#            if Imax > 0:
3263#                WriteCIFitem(self.fp, '   _gsas_i100_meas')
3264
3265            refcount = 0
3266            hklmin = None
3267            hklmax = None
3268            dmax = None
3269            dmin = None
3270            for phasenam in phaselist:
3271                scale = self.Phases[phasenam]['Histograms'][histlbl]['Scale'][0]
3272                phaseid = self.Phases[phasenam]['pId']
3273                refcount += len(histblk['Reflection Lists'][phasenam]['RefList'])
3274                refList = np.asarray(histblk['Reflection Lists'][phasenam]['RefList'])
3275                I100 = scale*refList.T[8]*refList.T[11]
3276                for j,ref in enumerate(histblk['Reflection Lists'][phasenam]['RefList']):
3277                    if DEBUG:
3278                        print('DEBUG: skipping reflection list')
3279                        break
3280                    if hklmin is None:
3281                        hklmin = copy.copy(ref[0:3])
3282                        hklmax = copy.copy(ref[0:3])
3283                    if dmin is None:
3284                         dmax = dmin = ref[4]
3285                    if len(phaselist) > 1:
3286                        s = PutInCol(phaseid,2)
3287                    else:
3288                        s = ""
3289                    for i,hkl in enumerate(ref[0:3]):
3290                        hklmax[i] = max(hkl,hklmax[i])
3291                        hklmin[i] = min(hkl,hklmin[i])
3292                        s += PutInCol(int(hkl),4)
3293                    for I in ref[8:10]:
3294                        s += PutInCol(G2mth.ValEsd(I,-0.0009),14)
3295                    s += PutInCol(G2mth.ValEsd(ref[10],-0.9),7)
3296                    dmax = max(dmax,ref[4])
3297                    dmin = min(dmin,ref[4])
3298                    s += PutInCol(G2mth.ValEsd(ref[4],-0.00009),8)
3299#                    if Imax > 0:
3300#                        s += PutInCol(G2mth.ValEsd(100.*I100[j]/Imax,-0.09),6)
3301                    s += PutInCol('o',2)
3302                    s += PutInCol('1',2)
3303                    s += PutInCol('1',2)
3304                    s += PutInCol('1',2)
3305                    s += PutInCol('.',2)
3306                    WriteCIFitem(self.fp, "  "+s)
3307
3308            # Write reflection statistics
3309            WriteCIFitem(self.fp, '_diffrn_reflns.number', str(refcount))
3310            if hklmin is not None and len(phaselist) == 1: # hkl range has no meaning with multiple phases
3311                WriteCIFitem(self.fp, '_diffrn_reflns.limit_h_min', str(int(hklmin[0])))
3312                WriteCIFitem(self.fp, '_diffrn_reflns.limit_h_max', str(int(hklmax[0])))
3313                WriteCIFitem(self.fp, '_diffrn_reflns.limit_k_min', str(int(hklmin[1])))
3314                WriteCIFitem(self.fp, '_diffrn_reflns.limit_k_max', str(int(hklmax[1])))
3315                WriteCIFitem(self.fp, '_diffrn_reflns.limit_l_min', str(int(hklmin[2])))
3316                WriteCIFitem(self.fp, '_diffrn_reflns.limit_l_max', str(int(hklmax[2])))
3317            if hklmin is not None:
3318                WriteCIFitem(self.fp, '_reflns.d_resolution_low  ', G2mth.ValEsd(dmax,-0.009))
3319                WriteCIFitem(self.fp, '_reflns.d_resolution_high ', G2mth.ValEsd(dmin,-0.009))
3320               
3321            WriteCIFitem(self.fp, '\n# POWDER DATA TABLE')
3322
3323            # is data fixed step? If the step varies by <0.01% treat as fixed step
3324            fixedstep = False
3325            zero = None
3326            WriteCIFitem(self.fp, '_refine.pdbx_pd_meas_number_of_points', str(len(histblk['Data'][0])))
3327            WriteCIFitem(self.fp, '\nloop_')
3328            #            WriteCIFitem(self.fp, '   _pd_proc_d_spacing') # need easy way to get this
3329            if 'T' in inst['Type'][0]: # and not TOF
3330                WriteCIFitem(self.fp, '   _pd_meas_time_of_flight')
3331            else:
3332                WriteCIFitem(self.fp, '   _pdbx_powder_data.pd_meas_2theta_scan')
3333            # at least for now, always report weights.
3334            #if countsdata:
3335            #    WriteCIFitem(self.fp, '   _pd_meas_counts_total')
3336            #else:
3337            WriteCIFitem(self.fp, '   _pdbx_powder_data.pd_meas_intensity_total')
3338            WriteCIFitem(self.fp, '   _pdbx_powder_data.pd_calc_intensity_total')
3339            WriteCIFitem(self.fp, '   _pdbx_powder_data.pd_proc_intensity_bkg_calc')
3340            WriteCIFitem(self.fp, '   _pdbx_powder_data.pd_proc_ls_weight')
3341            maxY = max(histblk['Data'][1].max(),histblk['Data'][3].max())
3342            if maxY < 0: maxY *= -10 # this should never happen, but...
3343            ndec = max(0,10-int(np.log10(maxY))-1) # 10 sig figs should be enough
3344            maxSU = histblk['Data'][2].max()
3345            if maxSU < 0: maxSU *= -1 # this should never happen, but...
3346            ndecSU = max(0,8-int(np.log10(maxSU))-1) # 8 sig figs should be enough
3347            lowlim,highlim = histblk['Limits'][1]
3348
3349            if DEBUG:
3350                print('DEBUG: skipping profile list')
3351            else:
3352                for x,yobs,yw,ycalc,ybkg in zip(histblk['Data'][0].data,        #get the data from these masked arrays
3353                                                histblk['Data'][1].data,
3354                                                histblk['Data'][2].data,
3355                                                histblk['Data'][3].data,
3356                                                histblk['Data'][4].data):
3357                    if lowlim <= x <= highlim:
3358                        pass
3359                    else:
3360                        yw = 0.0 # show the point is not in use
3361
3362                    if fixedstep:
3363                        s = ""
3364                    elif zero:
3365                        s = PutInCol(G2mth.ValEsd(x-zero,-0.00009),10)
3366                    else:
3367                        s = PutInCol(G2mth.ValEsd(x,-0.00009),10)
3368                    s += PutInCol(Yfmt(ndec,yobs),12)
3369                    s += PutInCol(Yfmt(ndec,ycalc),12)
3370                    s += PutInCol(Yfmt(ndec,ybkg),11)
3371                    s += PutInCol(Yfmt(ndecSU,yw),9)
3372                    WriteCIFitem(self.fp, "  "+s)
3373
3374        def WriteSingleXtalData(histlbl):
3375            'Write out the selected single crystal histogram info'
3376            histblk = self.Histograms[histlbl]
3377
3378            #refprx = '_refln.' # mm
3379            refprx = '_refln_' # normal
3380
3381            WriteCIFitem(self.fp, '\n# STRUCTURE FACTOR TABLE')
3382            WriteCIFitem(self.fp, 'loop_' +
3383                         '\n   ' + refprx + 'index_h' +
3384                         '\n   ' + refprx + 'index_k' +
3385                         '\n   ' + refprx + 'index_l' +
3386                         '\n   ' + refprx + 'F_squared_meas' +
3387                         '\n   ' + refprx + 'F_squared_sigma' +
3388                         '\n   ' + refprx + 'F_squared_calc' +
3389                         '\n   ' + refprx + 'phase_calc'
3390                         )
3391
3392            hklmin = None
3393            hklmax = None
3394            dmax = None
3395            dmin = None
3396            refcount = len(histblk['Data']['RefList'])
3397            for ref in histblk['Data']['RefList']:
3398                if ref[3] <= 0:      #skip user rejected reflections (mul <= 0)
3399                    continue
3400                s = "  "
3401                if hklmin is None:
3402                    hklmin = copy.copy(ref[0:3])
3403                    hklmax = copy.copy(ref[0:3])
3404                    dmax = dmin = ref[4]
3405                for i,hkl in enumerate(ref[0:3]):
3406                    hklmax[i] = max(hkl,hklmax[i])
3407                    hklmin[i] = min(hkl,hklmin[i])
3408                    s += PutInCol(int(hkl),4)
3409                if ref[5] == 0.0:
3410                    s += PutInCol(G2mth.ValEsd(ref[8],0),12)
3411                    s += PutInCol('.',10)
3412                    s += PutInCol(G2mth.ValEsd(ref[9],0),12)
3413                else:
3414                    sig = ref[6] * ref[8] / ref[5]
3415                    s += PutInCol(G2mth.ValEsd(ref[8],-abs(sig/10)),12)
3416                    s += PutInCol(G2mth.ValEsd(sig,-abs(sig)/10.),10)
3417                    s += PutInCol(G2mth.ValEsd(ref[9],-abs(sig/10)),12)
3418                s += PutInCol(G2mth.ValEsd(ref[10],-0.9),7)
3419                dmax = max(dmax,ref[4])
3420                dmin = min(dmin,ref[4])
3421                WriteCIFitem(self.fp, s)
3422            if not self.quickmode: # statistics only in a full CIF
3423                WriteReflStat(refcount,hklmin,hklmax,dmin,dmax)
3424                hId = histblk['hId']
3425                hfx = '0:'+str(hId)+':'
3426                phfx = '%d:%d:'%(0,hId)
3427                extType,extModel,extParms = self.Phases[phasenam]['Histograms'][histlbl]['Extinction']
3428                if extModel != 'None':
3429                    WriteCIFitem(self.fp, '# Extinction scaled by 1.e5')
3430                    WriteCIFitem(self.fp, '_refine_ls_extinction_method','Becker-Coppens %s %s'%(extModel,extType))
3431                    sig = -1.e-3
3432                    if extModel == 'Primary':
3433                        parm = extParms['Ep'][0]*1.e5
3434                        if extParms['Ep'][1]:
3435                            sig = self.sigDict[phfx+'Ep']*1.e5
3436                        text = G2mth.ValEsd(parm,sig)
3437                    elif extModel == 'Secondary Type I':
3438                        parm = extParms['Eg'][0]*1.e5
3439                        if extParms['Eg'][1]:
3440                            sig = self.sigDict[phfx+'Eg']*1.e5
3441                        text = G2mth.ValEsd(parm,sig)
3442                    elif extModel == 'Secondary Type II':
3443                        parm = extParms['Es'][0]*1.e5
3444                        if extParms['Es'][1]:
3445                            sig = self.sigDict[phfx+'Es']*1.e5
3446                        text = G2mth.ValEsd(parm,sig)
3447                    elif extModel == 'Secondary Type I & II':
3448                        parm = extParms['Eg'][0]*1.e5
3449                        if extParms['Es'][1]:
3450                            sig = self.sigDict[phfx+'Es']*1.e5
3451                        text = G2mth.ValEsd(parm,sig)
3452                        sig = -1.0e-3
3453                        parm = extParms['Es'][0]*1.e5
3454                        if extParms['Es'][1]:
3455                            sig = self.sigDict[phfx+'Es']*1.e5
3456                        text += G2mth.ValEsd(parm,sig)
3457                    WriteCIFitem(self.fp, '_refine_ls_extinction_coef',text)
3458                    WriteCIFitem(self.fp, '_refine_ls_extinction_expression','Becker & Coppens (1974). Acta Cryst. A30, 129-147')
3459
3460                WriteCIFitem(self.fp, '_refine_ls_wR_factor_gt    ','%.4f'%(histblk['wR']/100.))
3461                WriteCIFitem(self.fp, '_refine_ls_R_factor_gt     ','%.4f'%(histblk[hfx+'Rf']/100.))
3462                WriteCIFitem(self.fp, '_refine_ls_R_Fsqd_factor   ','%.4f'%(histblk[hfx+'Rf^2']/100.))
3463        def EditAuthor(event=None):
3464            'dialog to edit the CIF author info'
3465            'Edit the CIF author name'
3466            dlg = G2G.SingleStringDialog(self.G2frame,
3467                                          'Get CIF Author',
3468                                          'Provide CIF Author name (Last, First)',
3469                                          value=self.author)
3470            if not dlg.Show():
3471                dlg.Destroy()
3472                return False  # cancel was pressed
3473            self.author = dlg.GetValue()
3474            self.shortauthorname = self.author.replace(',','').replace(' ','')[:20]
3475            dlg.Destroy()
3476            try:
3477                self.OverallParms['Controls']["Author"] = self.author # save for future
3478            except KeyError:
3479                pass
3480            return True
3481
3482        def EditInstNames(event=None):
3483            'Provide a dialog for editing instrument names; for sequential fit, only need one name'
3484            dictlist = []
3485            keylist = []
3486            lbllist = []
3487            for hist in sorted(self.Histograms):
3488                if hist.startswith("PWDR"):
3489                    key2 = "Sample Parameters"
3490                    d = self.Histograms[hist][key2]
3491                elif hist.startswith("HKLF"):
3492                    key2 = "Instrument Parameters"
3493                    d = self.Histograms[hist][key2][0]
3494
3495                lbllist.append(hist)
3496                dictlist.append(d)
3497                keylist.append('InstrName')
3498                instrname = d.get('InstrName')
3499                if instrname is None:
3500                    d['InstrName'] = ''
3501                if hist.startswith("PWDR") and seqmode: break               
3502            return G2G.CallScrolledMultiEditor(
3503                self.G2frame,dictlist,keylist,
3504                prelbl=range(1,len(dictlist)+1),
3505                postlbl=lbllist,
3506                title='Instrument names',
3507                header="Edit instrument names. Note that a non-blank\nname is required for all histograms",
3508                CopyButton=True,ASCIIonly=True)
3509
3510        def EditRanges(event):
3511            '''Edit the bond distance/angle search range; phase is determined from
3512            a pointer placed in the button object (.phasedict) that references the
3513            phase dictionary
3514            '''
3515            but = event.GetEventObject()
3516            phasedict = but.phasedict
3517            dlg = G2G.DisAglDialog(
3518                self.G2frame,
3519                phasedict['General']['DisAglCtls'], # edited
3520                phasedict['General'], # defaults
3521                )
3522            if dlg.ShowModal() == wx.ID_OK:
3523                phasedict['General']['DisAglCtls'] = dlg.GetData()
3524            dlg.Destroy()
3525           
3526        def SetCellT(event):
3527            '''Set the temperature value by selection of a histogram
3528            '''
3529            but = event.GetEventObject()
3530            phasenam = but.phase
3531            rId =  self.Phases[phasenam]['ranId']
3532            self.CellHistSelection[rId] = self._CellSelectHist(phasenam)
3533           
3534        def EditCIFDefaults():
3535            '''Fills the CIF Defaults window with controls for editing various CIF export
3536            parameters (mostly related to templates).
3537            '''
3538            if len(self.cifdefs.GetChildren()) > 0:
3539                saveSize = self.cifdefs.GetSize()
3540                self.cifdefs.DestroyChildren()
3541            else:
3542                saveSize = None
3543            self.cifdefs.SetTitle('Edit CIF settings')
3544            vbox = wx.BoxSizer(wx.VERTICAL)
3545            vbox.Add(wx.StaticText(self.cifdefs, wx.ID_ANY,'Creating file '+self.filename))
3546            but = wx.Button(self.cifdefs, wx.ID_ANY,'Edit CIF Author')
3547            but.Bind(wx.EVT_BUTTON,EditAuthor)
3548            vbox.Add(but,0,wx.ALIGN_CENTER,3)
3549            but = wx.Button(self.cifdefs, wx.ID_ANY,'Edit Instrument Name(s)')
3550            but.Bind(wx.EVT_BUTTON,EditInstNames)
3551            vbox.Add(but,0,wx.ALIGN_CENTER,3)
3552            cpnl = wxscroll.ScrolledPanel(self.cifdefs,size=(300,300))
3553            cbox = wx.BoxSizer(wx.VERTICAL)
3554            G2G.HorizontalLine(cbox,cpnl)
3555            cbox.Add(
3556                CIFtemplateSelect(self.cifdefs,
3557                                  cpnl,'publ',self.OverallParms['Controls'],
3558                                  EditCIFDefaults,
3559                                  "Publication (overall) template",
3560                                  ),
3561                0,wx.EXPAND|wx.ALIGN_LEFT|wx.ALL)
3562            for phasenam in sorted(self.Phases.keys()):
3563                G2G.HorizontalLine(cbox,cpnl)
3564                title = 'Phase '+phasenam
3565                phasedict = self.Phases[phasenam] # pointer to current phase info
3566                cbox.Add(
3567                    CIFtemplateSelect(self.cifdefs,
3568                                      cpnl,'phase',phasedict['General'],
3569                                      EditCIFDefaults,
3570                                      title,
3571                                      phasenam),
3572                    0,wx.EXPAND|wx.ALIGN_LEFT|wx.ALL)
3573                cpnl.SetSizer(cbox)
3574                if phasedict['General']['Type'] == 'nuclear':
3575                    but = wx.Button(cpnl, wx.ID_ANY,'Edit distance/angle ranges')
3576                    cbox.Add(but,0,wx.ALIGN_LEFT,0)
3577                    cbox.Add((-1,2))
3578                    but.phasedict = self.Phases[phasenam]  # set a pointer to current phase info
3579                    but.Bind(wx.EVT_BUTTON,EditRanges)     # phase bond/angle ranges
3580                    but = wx.Button(cpnl, wx.ID_ANY,'Set distance/angle publication flags')
3581                    but.phase = phasenam  # set a pointer to current phase info
3582                    but.Bind(wx.EVT_BUTTON,SelectDisAglFlags)     # phase bond/angle ranges
3583                    cbox.Add(but,0,wx.ALIGN_LEFT,0)
3584                if self._CellSelectNeeded(phasenam):
3585                    but = wx.Button(cpnl, wx.ID_ANY,'Select cell temperature')
3586                    cbox.Add(but,0,wx.ALIGN_LEFT,0)
3587                    cbox.Add((-1,2))
3588                    but.phase = phasenam  # set a pointer to current phase info
3589                    but.Bind(wx.EVT_BUTTON,SetCellT)
3590                cbox.Add((-1,2))
3591            for i in sorted(self.powderDict.keys()):
3592                G2G.HorizontalLine(cbox,cpnl)
3593                if seqmode:
3594                    hist = self.powderDict[i]
3595                    histblk = self.Histograms[hist]
3596                    title = 'All Powder datasets'
3597                    cbox.Add(
3598                        CIFtemplateSelect(self.cifdefs,
3599                                      cpnl,'powder',self.OverallParms['Controls'],
3600                                      EditCIFDefaults,
3601                                      title,
3602                                      histblk["Sample Parameters"]['InstrName'],
3603                                      cifKey="seqCIF_template"),
3604                        0,wx.EXPAND|wx.ALIGN_LEFT|wx.ALL)
3605                    break
3606                hist = self.powderDict[i]
3607                histblk = self.Histograms[hist]
3608                title = 'Powder dataset '+hist[5:]
3609                cbox.Add(
3610                    CIFtemplateSelect(self.cifdefs,
3611                                      cpnl,'powder',histblk["Sample Parameters"],
3612                                      EditCIFDefaults,
3613                                      title,
3614                                      histblk["Sample Parameters"]['InstrName']),
3615                    0,wx.EXPAND|wx.ALIGN_LEFT|wx.ALL)
3616            for i in sorted(self.xtalDict.keys()):
3617                G2G.HorizontalLine(cbox,cpnl)
3618                hist = self.xtalDict[i]
3619                histblk = self.Histograms[hist]
3620                title = 'Single Xtal dataset '+hist[5:]
3621                cbox.Add(
3622                    CIFtemplateSelect(self.cifdefs,
3623                                      cpnl,'single',histblk["Instrument Parameters"][0],
3624                                      EditCIFDefaults,
3625                                      title,
3626                                      histblk["Instrument Parameters"][0]['InstrName']),
3627                    0,wx.EXPAND|wx.ALIGN_LEFT|wx.ALL)
3628            cpnl.SetSizer(cbox)
3629            cpnl.SetAutoLayout(1)
3630            cpnl.SetupScrolling()
3631            #cpnl.Bind(rw.EVT_RW_LAYOUT_NEEDED, self.OnLayoutNeeded) # needed if sizes change
3632            cpnl.Layout()
3633
3634            vbox.Add(cpnl, 1, wx.ALIGN_LEFT|wx.ALL|wx.EXPAND, 0)
3635            btnsizer = wx.StdDialogButtonSizer()
3636            btn = wx.Button(self.cifdefs, wx.ID_OK, "Create CIF")
3637            btn.SetDefault()
3638            btnsizer.AddButton(btn)
3639            btn = wx.Button(self.cifdefs, wx.ID_CANCEL)
3640            btnsizer.AddButton(btn)
3641            btnsizer.Realize()
3642            vbox.Add(btnsizer, 0, wx.ALIGN_CENTER|wx.ALL, 5)
3643            self.cifdefs.SetSizer(vbox)
3644            if not saveSize:
3645                vbox.Fit(self.cifdefs)
3646            self.cifdefs.Layout()
3647
3648        def OnToggleButton(event):
3649            'Respond to press of ToggleButton in SelectDisAglFlags'
3650            but = event.GetEventObject()
3651            if but.GetValue():
3652                but.DisAglSel[but.key] = True
3653            else:
3654                try:
3655                    del but.DisAglSel[but.key]
3656                except KeyError:
3657                    pass
3658        def keepTrue(event):
3659            event.GetEventObject().SetValue(True)
3660        def keepFalse(event):
3661            event.GetEventObject().SetValue(False)
3662
3663        def SelectDisAglFlags(event):
3664            'Select Distance/Angle use flags for the selected phase'
3665            phasenam = event.GetEventObject().phase
3666            phasedict = self.Phases[phasenam]
3667            SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(phasedict['General']['SGData'])
3668            generalData = phasedict['General']
3669            # create a dict for storing Pub flag for bonds/angles, if needed
3670            if phasedict['General'].get("DisAglHideFlag") is None:
3671                phasedict['General']["DisAglHideFlag"] = {}
3672            DisAngSel = phasedict['General']["DisAglHideFlag"]
3673
3674            cx,ct,cs,cia = phasedict['General']['AtomPtrs']
3675            cn = ct-1
3676            cfrac = cx+3
3677            DisAglData = {}
3678            # create a list of atoms, but skip atoms with zero occupancy
3679            xyz = []
3680            fpfx = str(phasedict['pId'])+'::Afrac:'
3681            for i,atom in enumerate(phasedict['Atoms']):
3682                if self.parmDict.get(fpfx+str(i),atom[cfrac]) == 0.0: continue
3683                xyz.append([i,]+atom[cn:cn+2]+atom[cx:cx+3])
3684            if 'DisAglCtls' not in generalData:
3685                # should not be used, since DisAglDialog should be called
3686                # for all phases before getting here
3687                dlg = G2G.DisAglDialog(
3688                    self.cifdefs,
3689                    {},
3690                    generalData)
3691                if dlg.ShowModal() == wx.ID_OK:
3692                    generalData['DisAglCtls'] = dlg.GetData()
3693                else:
3694                    dlg.Destroy()
3695                    return
3696                dlg.Destroy()
3697            dlg = wx.Dialog(
3698                self.G2frame,
3699                style=wx.DEFAULT_DIALOG_STYLE | wx.RESIZE_BORDER)
3700            vbox = wx.BoxSizer(wx.VERTICAL)
3701            txt = wx.StaticText(dlg,wx.ID_ANY,'Searching distances for phase '+phasenam
3702                                +'\nPlease wait...')
3703            vbox.Add(txt,0,wx.ALL|wx.EXPAND)
3704            dlg.SetSizer(vbox)
3705            dlg.CenterOnParent()
3706            dlg.Show() # post "please wait"
3707            wx.BeginBusyCursor() # and change cursor
3708
3709            DisAglData['OrigAtoms'] = xyz
3710            DisAglData['TargAtoms'] = xyz
3711            SymOpList,offsetList,symOpList,G2oprList,G2opcodes = G2spc.AllOps(
3712                generalData['SGData'])
3713
3714#            xpandSGdata = generalData['SGData'].copy()
3715#            xpandSGdata.update({'SGOps':symOpList,
3716#                                'SGInv':False,
3717#                                'SGLatt':'P',
3718#                                'SGCen':np.array([[0, 0, 0]]),})
3719#            DisAglData['SGData'] = xpandSGdata
3720            DisAglData['SGData'] = generalData['SGData'].copy()
3721
3722            DisAglData['Cell'] = generalData['Cell'][1:] #+ volume
3723            if 'pId' in phasedict:
3724                DisAglData['pId'] = phasedict['pId']
3725                DisAglData['covData'] = self.OverallParms['Covariance']
3726            try:
3727                AtomLabels,DistArray,AngArray = G2stMn.RetDistAngle(
3728                    generalData['DisAglCtls'],
3729                    DisAglData)
3730            except KeyError:        # inside DistAngle for missing atom types in DisAglCtls
3731                print('**** ERROR - try again but do "Reset" to fill in missing atom types ****')
3732            wx.EndBusyCursor()
3733            txt.SetLabel('Set publication flags for distances and angles in\nphase '+phasenam)
3734            vbox.Add((5,5))
3735            vbox.Add(wx.StaticText(dlg,wx.ID_ANY,
3736                                   'The default is to flag all distances and angles as to be'+
3737                                   '\npublished. Change this by pressing appropriate buttons.'),
3738                     0,wx.ALL|wx.EXPAND)
3739            hbox = wx.BoxSizer(wx.HORIZONTAL)
3740            vbox.Add(hbox)
3741            hbox.Add(wx.StaticText(dlg,wx.ID_ANY,'Button appearance: '))
3742            but = wx.ToggleButton(dlg,wx.ID_ANY,'Publish')
3743            but.Bind(wx.EVT_TOGGLEBUTTON,keepFalse)
3744            hbox.Add(but)
3745            but = wx.ToggleButton(dlg,wx.ID_ANY,"Don't publish")
3746            but.Bind(wx.EVT_TOGGLEBUTTON,keepTrue)
3747            hbox.Add(but)
3748            but.SetValue(True)
3749            G2G.HorizontalLine(vbox,dlg)
3750
3751            cpnl = wxscroll.ScrolledPanel(dlg,size=(400,300))
3752            cbox = wx.BoxSizer(wx.VERTICAL)
3753            for c in sorted(DistArray):
3754                karr = []
3755                UsedCols = {}
3756                cbox.Add(wx.StaticText(cpnl,wx.ID_ANY,
3757                                   'distances to/angles around atom '+AtomLabels[c]))
3758                #dbox = wx.GridBagSizer(hgap=5)
3759                dbox = wx.GridBagSizer()
3760                for i,D in enumerate(DistArray[c]):
3761                    karr.append(tuple(D[0:3]))
3762                    val = "{:.2f}".format(D[3])
3763                    sym = " [{:d} {:d} {:d}]".format(*D[1]) + " #{:d}".format(D[2])
3764                    dbox.Add(wx.StaticText(cpnl,wx.ID_ANY,AtomLabels[D[0]]),
3765                             (i+1,0)
3766                             )
3767                    dbox.Add(wx.StaticText(cpnl,wx.ID_ANY,sym),
3768                             (i+1,1)
3769                             )
3770                    but = wx.ToggleButton(cpnl,wx.ID_ANY,val)
3771                    but.key = (c,karr[-1])
3772                    but.DisAglSel = DisAngSel
3773                    if DisAngSel.get(but.key): but.SetValue(True)
3774                    but.Bind(wx.EVT_TOGGLEBUTTON,OnToggleButton)
3775                    dbox.Add(but,(i+1,2),border=1)
3776                for i,D in enumerate(AngArray[c]):
3777                    val = "{:.1f}".format(D[2][0])
3778                    but = wx.ToggleButton(cpnl,wx.ID_ANY,val)
3779                    but.key = (karr[D[0]],c,karr[D[1]])
3780                    but.DisAglSel = DisAngSel
3781                    if DisAngSel.get(but.key): but.SetValue(True)
3782                    but.Bind(wx.EVT_TOGGLEBUTTON,OnToggleButton)
3783                    dbox.Add(but,(D[0]+1,D[1]+3),border=1)
3784                    UsedCols[D[1]+3] = True
3785                for i,D in enumerate(DistArray[c][:-1]): # label columns that are used
3786                    if UsedCols.get(i+3):
3787                        dbox.Add(wx.StaticText(cpnl,wx.ID_ANY,AtomLabels[D[0]]),
3788                                 (0,i+3),
3789                                 flag=wx.ALIGN_CENTER
3790                                 )
3791                dbox.Add(wx.StaticText(cpnl,wx.ID_ANY,'distance'),
3792                                 (0,2),
3793                                 flag=wx.ALIGN_CENTER
3794                                 )
3795                cbox.Add(dbox)
3796                G2G.HorizontalLine(cbox,cpnl)
3797            cpnl.SetSizer(cbox)
3798            cpnl.SetAutoLayout(1)
3799            cpnl.SetupScrolling()
3800            #cpnl.Bind(rw.EVT_RW_LAYOUT_NEEDED, self.OnLayoutNeeded) # needed if sizes change
3801            cpnl.Layout()
3802
3803            vbox.Add(cpnl, 1, wx.ALIGN_LEFT|wx.ALL|wx.EXPAND, 0)
3804
3805            btnsizer = wx.StdDialogButtonSizer()
3806            btn = wx.Button(dlg, wx.ID_OK, "Done")
3807            btn.SetDefault()
3808            btnsizer.AddButton(btn)
3809            btnsizer.Realize()
3810            vbox.Add(btnsizer, 0, wx.ALIGN_CENTER|wx.ALL, 5)
3811            dlg.SetSizer(vbox)
3812            vbox.Fit(dlg)
3813            dlg.Layout()
3814
3815            dlg.CenterOnParent()
3816            dlg.ShowModal()
3817
3818#==============================================================================
3819####  _Exporter code starts here         ======================================
3820#==============================================================================
3821        # make sure required information is present
3822        self.CIFdate = dt.datetime.strftime(dt.datetime.now(),"%Y-%m-%dT%H:%M")
3823        if not self.CIFname: # Get a name for the CIF. If not defined, use the GPX name (save, if that is needed).
3824            if not self.G2frame.GSASprojectfile:
3825                self.G2frame.OnFileSaveas(None)
3826            if not self.G2frame.GSASprojectfile: return
3827            self.CIFname = os.path.splitext(
3828                os.path.split(self.G2frame.GSASprojectfile)[1]
3829                )[0]
3830            self.CIFname = self.CIFname.replace(' ','')
3831        # replace non-ASCII characters in CIFname with dots
3832        s = ''
3833        for c in self.CIFname:
3834            if ord(c) < 128:
3835                s += c
3836            else:
3837                s += '.'
3838        self.CIFname = s
3839        phasebyhistDict = {} # a cross-reference to phases by histogram -- sequential fits
3840        #=================================================================
3841        # write quick CIFs
3842        #=================================================================
3843        if phaseOnly: #====Phase only CIF ================================
3844            print('Writing CIF output to file '+self.filename)
3845            oneblock = True
3846            self.quickmode = True
3847            self.Write(' ')
3848            self.Write(70*'#')
3849            WriteCIFitem(self.fp, 'data_'+phaseOnly.replace(' ','_'))
3850            #phaseblk = self.Phases[phaseOnly] # pointer to current phase info
3851            # report the phase info
3852            if self.Phases[phaseOnly]['General']['Type'] == 'macromolecular':
3853                WritePhaseInfoMM(phaseOnly)
3854            else:
3855                WritePhaseInfo(phaseOnly)
3856            return
3857        elif histOnly: #====Histogram only CIF ================================
3858            print('Writing CIF output to file '+self.filename)
3859            MM = False
3860            for p in self.Phases:
3861                if self.Phases[p]['General']['Type'] == 'macromolecular':
3862                    MM = True
3863                    break
3864            hist = histOnly
3865            #histname = histOnly.replace(' ','')
3866            oneblock = True
3867            self.quickmode = True
3868            self.ifHKLF = False
3869            self.ifPWDR = True
3870            self.Write(' ')
3871            self.Write(70*'#')
3872            #phasenam = self.Phases.keys()[0]
3873            WriteCIFitem(self.fp, 'data_'+self.CIFname)
3874            if hist.startswith("PWDR") and MM:
3875                WritePowderDataMM(hist)
3876            elif hist.startswith("PWDR"):
3877                WritePowderData(hist)
3878            elif hist.startswith("HKLF"):
3879                WriteSingleXtalData(hist)
3880            return
3881        #===============================================================================
3882        # setup for sequential fits here
3883        #===============================================================================
3884        seqmode = False
3885        seqHistList = []
3886        if self.G2frame.testSeqRefineMode():
3887            if self.seqData is None:
3888                raise Exception('Use Export/Sequential project for sequential refinements')
3889            if len(self.Phases) > 1:
3890                phaseWithHist = False  # multiple phases per histogram
3891            else:
3892                phaseWithHist = True   # include the phase in the same block as the histogram
3893            seqmode = True
3894            seqHistList = [h for h in self.seqData['histNames'] if h in self.seqData]
3895            if 'Use' in self.seqData and len(seqHistList) == len(self.seqData.get('Use',[])):
3896                seqHistList = [h for i,h in enumerate(seqHistList) if self.seqData['Use'][i]]
3897               
3898        #===============================================================================
3899        ### full CIF export starts here (Sequential too)
3900        #===============================================================================
3901        # load saved CIF author name
3902        self.author = self.OverallParms['Controls'].get("Author",'?').strip()
3903        # initialize dict for Selection of Hist for unit cell reporting
3904        self.OverallParms['Controls']['CellHistSelection'] = self.OverallParms[
3905            'Controls'].get('CellHistSelection',{})
3906        self.CellHistSelection = self.OverallParms['Controls']['CellHistSelection']
3907        # create a dict with refined values and their uncertainties
3908        self.loadParmDict()
3909        # is there anything to export?
3910        if len(self.Phases) == len(self.powderDict) == len(self.xtalDict) == 0:
3911           self.G2frame.ErrorDialog(
3912               'Empty project',
3913               'Project does not contain any data or phases. Are they interconnected?')
3914           return
3915        if self.ExportSelect('ask'): return
3916        if not self.filename:
3917            print('No name supplied')
3918            return
3919        self.OpenFile(delayOpen=True)
3920       
3921        #if self.ExportSelect('default'): return
3922        # Someday: get restraint & constraint info
3923        #restraintDict = self.OverallParms.get('Restraints',{})
3924        #for i in  self.OverallParms['Constraints']:
3925        #    print i
3926        #    for j in self.OverallParms['Constraints'][i]:
3927        #        print j
3928
3929        self.quickmode = False # full CIF
3930        phasenam = None # include all phases
3931        # Will this require a multiblock CIF?
3932        if len(self.Phases) > 1:
3933            oneblock = False
3934        elif len(self.powderDict) + len(self.xtalDict) > 1:
3935            oneblock = False
3936        else: # one phase, one dataset, Full CIF
3937            oneblock = True
3938
3939        # check there is an instrument name for every histogram
3940        self.ifPWDR = False
3941        self.ifHKLF = False
3942        invalid = 0
3943        key3 = 'InstrName'
3944        for hist in self.Histograms:
3945            if hist.startswith("PWDR"):
3946                self.ifPWDR = True
3947                key2 = "Sample Parameters"
3948                d = self.Histograms[hist][key2]
3949            elif hist.startswith("HKLF"):
3950                self.ifHKLF = True
3951                key2 = "Instrument Parameters"
3952                d = self.Histograms[hist][key2][0]
3953            instrname = d.get(key3)
3954            if instrname is None:
3955                d[key3] = ''
3956                invalid += 1
3957            elif instrname.strip() == '':
3958                invalid += 1
3959            if hist.startswith("PWDR") and seqmode: break
3960        if invalid:
3961            #msg = ""
3962            #if invalid > 3: msg = (
3963            #    "\n\nNote: it may be faster to set the name for\n"
3964            #    "one histogram for each instrument and use the\n"
3965            #    "File/Copy option to duplicate the name"
3966            #    )
3967            if not EditInstNames(): return
3968
3969        # check for a distance-angle range search range for each phase
3970        for phasenam in sorted(self.Phases.keys()):
3971            #i = self.Phases[phasenam]['pId']
3972            phasedict = self.Phases[phasenam] # pointer to current phase info
3973            if 'DisAglCtls' not in phasedict['General']:
3974                dlg = G2G.DisAglDialog(
3975                    self.G2frame,
3976                    {},
3977                    phasedict['General'])
3978                if dlg.ShowModal() == wx.ID_OK:
3979                    phasedict['General']['DisAglCtls'] = dlg.GetData()
3980                else:
3981                    dlg.Destroy()
3982                    return
3983                dlg.Destroy()
3984                   
3985        # check if temperature values & pressure are defaulted
3986        default = 0
3987        for hist in self.Histograms:
3988            if hist.startswith("PWDR"):
3989                key2 = "Sample Parameters"
3990                T = self.Histograms[hist][key2].get('Temperature')
3991                if not T:
3992                    default += 1
3993                elif T == 300:
3994                    default += 1
3995                P = self.Histograms[hist][key2].get('Pressure')
3996                if not P:
3997                    default += 1
3998                elif P == 1:
3999                    default += 1
4000        if default > 0:
4001            dlg = wx.MessageDialog(
4002                self.G2frame,
4003                'Temperature/Pressure values appear to be defaulted for some powder histograms (See Sample Parameters for each PWDR tree entry). Do you want to use those values?',
4004                'Check T and P values',
4005                wx.OK|wx.CANCEL)
4006            ret = dlg.ShowModal()
4007            dlg.CenterOnParent()
4008            dlg.Destroy()
4009            if ret != wx.ID_OK: return
4010        if oneblock:
4011            # select a dataset to use (there should only be one set in one block,
4012            # but take whatever comes 1st)
4013            for hist in self.Histograms:
4014                histblk = self.Histograms[hist]
4015                if hist.startswith("PWDR"):
4016                    instnam = histblk["Sample Parameters"]['InstrName']
4017                    break # ignore all but 1st data histogram
4018                elif hist.startswith("HKLF"):
4019                    instnam = histblk["Instrument Parameters"][0]['InstrName']
4020                    break # ignore all but 1st data histogram
4021        # give the user a window to edit CIF contents
4022        if not self.author:
4023            self.author = self.OverallParms['Controls'].get("Author",'?').strip()
4024        if not self.author:
4025            if not EditAuthor(): return
4026        self.ValidateAscii([('Author name',self.author),]) # check for ASCII strings where needed, warn on problems
4027        self.shortauthorname = self.author.replace(',','').replace(' ','')[:20]
4028        self.cifdefs = wx.Dialog(
4029            self.G2frame,
4030            style=wx.DEFAULT_DIALOG_STYLE | wx.RESIZE_BORDER)
4031        self.cifdefs.G2frame = self.G2frame
4032        self.cifdefs.CenterOnParent()
4033        EditCIFDefaults()
4034        if self.cifdefs.ShowModal() != wx.ID_OK:
4035            self.cifdefs.Destroy()
4036            return
4037        while self.ValidateAscii([('Author name',self.author),
4038                                  ]): # validate a few things as ASCII
4039            if self.cifdefs.ShowModal() != wx.ID_OK:
4040                self.cifdefs.Destroy()
4041                return
4042        self.cifdefs.Destroy()
4043        MM = False
4044        for p in self.Phases:
4045            if self.Phases[p]['General']['Type'] == 'macromolecular':
4046                MM = True
4047                break
4048        #======================================================================
4049        # export different types of CIFs below
4050        #======================================================================
4051        print('Writing CIF output to file '+self.filename+"...")
4052        self.openDelayed()
4053        if self.currentExportType == 'single' or self.currentExportType == 'powder':
4054            #======================================================================
4055            #### Data only CIF (powder/xtal) ======================================
4056            #======================================================================
4057            hist = self.histnam[0]
4058            self.CIFname = hist[5:40].replace(' ','')
4059            WriteCIFitem(self.fp, 'data_'+self.CIFname)
4060            if hist.startswith("PWDR") and MM:
4061                WritePowderDataMM(hist)
4062            elif hist.startswith("PWDR"):
4063                WritePowderData(hist)
4064            elif hist.startswith("HKLF"):
4065                WriteSingleXtalData(hist)
4066            else:
4067                print ("should not happen")
4068        elif oneblock:
4069            #======================================================================
4070            #### Full (data & phase) single block CIF =============================
4071            #======================================================================
4072            WriteCIFitem(self.fp, 'data_'+self.CIFname)
4073            if phasenam is None: # if not already selected, select the first phase (should be one)
4074                phasenam = self.Phases.keys()[0]
4075            #print 'phasenam',phasenam
4076            #phaseblk = self.Phases[phasenam] # pointer to current phase info
4077            instnam = instnam.replace(' ','')
4078            WriteCIFitem(self.fp, '_pd_block_id',
4079                         str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
4080                         str(self.shortauthorname) + "|" + instnam)
4081            WriteAudit()
4082            writeCIFtemplate(self.OverallParms['Controls'],'publ') # overall (publication) template
4083            if MM:
4084                WriteOverallMM()
4085            else:
4086                WriteOverall()
4087            writeCIFtemplate(self.Phases[phasenam]['General'],'phase',phasenam) # write phase template
4088            # report the phase info
4089            if self.Phases[phasenam]['General']['Type'] == 'macromolecular':
4090                WritePhaseInfoMM(phasenam,False)
4091            else:
4092                WritePhaseInfo(phasenam,False)
4093            if hist.startswith("PWDR"):  # this is invoked for single-block CIFs
4094                # preferred orientation
4095                SH = FormatSH(phasenam)
4096                MD = FormatHAPpo(phasenam)
4097                if SH and MD:
4098                    WriteCIFitem(self.fp, '_pd_proc_ls_pref_orient_corr', SH + '\n' + MD)
4099                elif SH or MD:
4100                    WriteCIFitem(self.fp, '_pd_proc_ls_pref_orient_corr', SH + MD)
4101                else:
4102                    WriteCIFitem(self.fp, '_pd_proc_ls_pref_orient_corr', 'none')
4103                # report profile, since one-block: include both histogram and phase info (N.B. there is only 1 of each)
4104                WriteCIFitem(self.fp, '_pd_proc_ls_profile_function',
4105                    FormatInstProfile(histblk["Instrument Parameters"],histblk['hId'])
4106                    +'\n'+FormatPhaseProfile(phasenam))
4107
4108                histblk = self.Histograms[hist]["Sample Parameters"]
4109                writeCIFtemplate(histblk,'powder',histblk['InstrName']) # write powder template
4110                if hist.startswith("PWDR") and MM:
4111                    WritePowderDataMM(hist)
4112                else:
4113                    WritePowderData(hist)
4114            elif hist.startswith("HKLF"):
4115                histprm = self.Histograms[hist]["Instrument Parameters"][0]
4116                writeCIFtemplate(histprm,'single',histprm['InstrName']) # single crystal template
4117                WriteSingleXtalData(hist)
4118        elif seqHistList:
4119            #======================================================================
4120            #### sequential fit export (multiblock)
4121            #======================================================================
4122            for phasenam in sorted(self.Phases.keys()):
4123                rId = phasedict['ranId']
4124                if rId in self.CellHistSelection: continue
4125                self.CellHistSelection[rId] = self._CellSelectHist(phasenam)
4126            nsteps = 1 + len(self.Phases) + len(seqHistList)
4127            try:
4128                dlg = wx.ProgressDialog('CIF progress','starting',nsteps,parent=self.G2frame)
4129                dlg.CenterOnParent()
4130
4131                # publication info block
4132                step = 1
4133                dlg.Update(step,"Exporting overall section")
4134                WriteCIFitem(self.fp, '\ndata_'+self.CIFname+'_publ')
4135                WriteAudit()
4136                WriteCIFitem(self.fp, '_pd_block_id',
4137                             str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
4138                             str(self.shortauthorname) + "|Overall")
4139                writeCIFtemplate(self.OverallParms['Controls'],'publ') #insert the publication template
4140                # ``template_publ.cif`` or a modified version
4141               
4142                # overall info block
4143                WriteCIFitem(self.fp, 'data_'+str(self.CIFname)+'_overall')
4144                WriteOverall('seq')
4145                hist = seqHistList[0]
4146                instnam = self.Histograms[hist]["Sample Parameters"]['InstrName']
4147                writeCIFtemplate(self.OverallParms['Controls'],'powder',instnam,
4148                                     cifKey="seqCIF_template") # powder template for all histograms
4149                instnam = instnam.replace(' ','')
4150                #============================================================
4151                if phaseWithHist:
4152                    WriteCIFitem(self.fp, '# POINTERS TO HISTOGRAM BLOCKS (Phase in histogram block)')
4153                else:
4154                    WriteCIFitem(self.fp, '# POINTERS TO HISTOGRAM BLOCKS (Phases pointer in histogram block)') 
4155                datablockidDict = {} # save block names here
4156                # loop over data blocks
4157                WriteCIFitem(self.fp, 'loop_   _pd_block_diffractogram_id')
4158                for hist in seqHistList:
4159                    j = self.Histograms[hist]['hId']
4160                    datablockidDict[hist] = (str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
4161                                             str(self.shortauthorname) + "|" +
4162                                             instnam + "_hist_"+str(j))
4163                    WriteCIFitem(self.fp, '  '+datablockidDict[hist])
4164                # for i in sorted(self.xtalDict.keys()):
4165                #     hist = self.xtalDict[i]
4166                #     histblk = self.Histograms[hist]
4167                #     instnam = histblk["Instrument Parameters"][0]['InstrName']
4168                #     instnam = instnam.replace(' ','')
4169                #     i = histblk['hId']
4170                #     datablockidDict[hist] = (str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
4171                #                              str(self.shortauthorname) + "|" +
4172                #                              instnam + "_hist_"+str(i))
4173                #     WriteCIFitem(self.fp, loopprefix,datablockidDict[hist])
4174                # setup and show sequential results table
4175                tblLabels,tblValues,tblSigs,tblTypes = mkSeqResTable('cif',seqHistList,self.seqData,
4176                                                    self.Phases,self.Histograms,self.Controls)
4177                WriteCIFitem(self.fp, '\n# Sequential results table') # (in case anyone can make sense of it)
4178                WriteCIFitem(self.fp, 'loop_   _gsas_seq_results_col_num _gsas_seq_results_col_label')
4179                for i,lbl in enumerate(tblLabels):
4180                    s = PutInCol(str(i),5)
4181                    if ' ' in lbl:
4182                        s += '"' + lbl + '"'
4183                    else:
4184                        s += lbl
4185                    WriteCIFitem(self.fp,"  "+s)
4186                s = 'loop_ '
4187                linelength = 120
4188                for i in range(len(tblLabels)):
4189                    if len(s) > linelength:
4190                        WriteCIFitem(self.fp,s)
4191                        s = '  '
4192                    s += " _gsas_seq_results_val" + str(i)
4193                WriteCIFitem(self.fp,s)
4194               
4195                for r in range(len(tblValues[0])):
4196                    s = ''
4197                    for c in range(len(tblLabels)):
4198                        if len(s) > linelength:
4199                            WriteCIFitem(self.fp,s)
4200                            s = '  '
4201                        sig = None
4202                        if tblSigs[c] is not None:
4203                            sig = tblSigs[c][r]
4204                           
4205                        if tblValues[c][r] is None:
4206                            if tblTypes[c] == 'int':
4207                                wid = 5
4208                            elif tblTypes[c] == 'str':
4209                                wid = 10
4210                            else:
4211                                wid = 12                               
4212                            s += PutInCol('.',wid)
4213                        elif sig is None and ',' in tblTypes[c]:
4214                            s += PutInCol(
4215                                ('{{:{}.{}f}}'.format(*tblTypes[c].split(','))).format(tblValues[c][r]),12)
4216                        elif tblTypes[c] == 'int':
4217                            s += PutInCol(str(tblValues[c][r]),5)
4218                        elif tblTypes[c] == 'str':
4219                            s += PutInCol(str(tblValues[c][r]),10)
4220                        elif sig is None and ',' in tblTypes[c]:
4221                            s += PutInCol(
4222                                ('{{:{}.{}f}}'.format(*tblTypes[c].split(','))).format(tblValues[c][r]),12)
4223                        elif sig is None and tblTypes[c] == 'float':
4224                            s += PutInCol('{:.6g}'.format(tblValues[c][r]),12)
4225                        elif sig:
4226                            s += PutInCol(G2mth.ValEsd(tblValues[c][r],sig),12)
4227                        else:
4228                            s += PutInCol(str(tblValues[c][r]),15)
4229                    WriteCIFitem(self.fp,s+'\n')
4230
4231                # sample template info: a block for each phase in project
4232                histblk = self.Histograms[seqHistList[0]]
4233                if phaseWithHist:    # include sample info in overall block
4234                    step += 1
4235                    dlg.Update(step,"Exporting phase")
4236                    phasenam = list(self.Phases.keys())[0]
4237                    writeCIFtemplate(self.Phases[phasenam]['General'],'phase',phasenam) # write phase template
4238                    WriteSeqOverallPhaseInfo(phasenam,histblk)
4239                else:
4240                    for j,phasenam in enumerate(sorted(self.Phases.keys())):
4241                        pId = self.Phases[phasenam]['pId']
4242                        step += 1
4243                        dlg.Update(step,"Exporting phase {}".format(pId))
4244                        WriteCIFitem(self.fp, '\n#'+78*'=')
4245                        WriteCIFitem(self.fp, 'data_'+self.CIFname+"_overall_phase"+str(j)+'\n')
4246                        writeCIFtemplate(self.Phases[phasenam]['General'],'phase',phasenam) # write phase template
4247                        WriteSeqOverallPhaseInfo(phasenam,histblk)
4248
4249                # create a block for each histogram, include phase in block for one-phase refinements
4250                # or separate blocks for each phase & histogram if more than one phase
4251                for i,hist in enumerate(seqHistList):
4252                    print('processing hist #',i,'hId=',self.Histograms[hist]['hId'],hist)
4253                    hId = self.Histograms[hist]['hId']
4254                    step += 1
4255                    dlg.Update(step,"Exporting "+hist.strip())
4256                    histblk = self.Histograms[hist]
4257                    WriteCIFitem(self.fp, '# Information for histogram '+str(i)+': '+hist)
4258                    WriteCIFitem(self.fp, '\ndata_'+self.CIFname+"_pwd_"+str(i))
4259                    WriteCIFitem(self.fp, '_pd_block_id',datablockidDict[hist])
4260                    if not phaseWithHist:
4261                        WriteCIFitem(self.fp, '\n# POINTERS TO PHASE BLOCKS')
4262                        phaseBlockName = {}
4263
4264                        wtFrSum = 0.
4265                        for j,phasenam in enumerate(sorted(self.Phases.keys())):
4266                            if hist not in self.Phases[phasenam]['Histograms']: continue
4267                            if not self.Phases[phasenam]['Histograms'][hist]['Use']: continue
4268                            phFrac = self.Phases[phasenam]['Histograms'][hist]['Scale'][0]
4269                            phFracKey = str(self.Phases[phasenam]['pId'])+':'+str(hId)+':Scale'
4270                            phFrac = self.seqData[hist]['parmDict'].get(phFracKey,phFrac)
4271                            wtFrSum += phFrac * self.Phases[phasenam]['General']['Mass']
4272                        WriteCIFitem(self.fp, 'loop_ _pd_phase_id _pd_phase_block_id _pd_phase_mass_%')
4273                        for j,phasenam in enumerate(sorted(self.Phases.keys())):
4274                            pId = self.Phases[phasenam]['pId']
4275                            if hist not in self.Phases[phasenam]['Histograms']: continue
4276                            if not self.Phases[phasenam]['Histograms'][hist]['Use']: continue
4277                            if ' ' in phasenam:
4278                                s = PutInCol('"'+phasenam+'"',20)
4279                            else:
4280                                s = PutInCol(phasenam,20)
4281                            phaseBlockName[pId] = datablockidDict[hist]+'_p'+str(j+1)
4282                            phFrac = self.Phases[phasenam]['Histograms'][hist]['Scale'][0]
4283                            phFracKey = str(self.Phases[phasenam]['pId'])+':'+str(hId)+':Scale'
4284                            phFrac = self.seqData[hist]['parmDict'].get(phFracKey,phFrac)
4285                            wtFr = phFrac * self.Phases[phasenam]['General']['Mass'] / wtFrSum
4286                            if phFracKey in self.seqData[hist]['varyList']:
4287                                sig = self.seqData[hist]['sig'][self.seqData[hist]['varyList'].index(phFracKey)]
4288                                sig *= self.Phases[phasenam]['General']['Mass'] / wtFrSum
4289                            elif phFracKey in self.seqData[hist]['depParmDict']:
4290                                sig = self.seqData[hist]['depParmDict'][phFracKey][1]
4291                                sig *= self.Phases[phasenam]['General']['Mass'] / wtFrSum
4292                            else:
4293                                sig = -0.0001
4294                            WriteCIFitem(self.fp, "  "+ s + " " + phaseBlockName[pId] + "  " + G2mth.ValEsd(wtFr,sig))
4295                        PP = FormatInstProfile(histblk["Instrument Parameters"],histblk['hId'])
4296                        PP += '\n'
4297                        WriteCIFitem(self.fp, '_pd_proc_ls_profile_function',PP)
4298               
4299                    WritePowderData(hist,seq=True) # write background, data & reflections, some instrument & sample terms
4300                    writeCIFtemplate(self.OverallParms['Controls'],'powder',
4301                                         self.Histograms[hist]["Sample Parameters"]['InstrName'],
4302                                         cifKey="seqCIF_template") # powder template for all histograms
4303                    WriteCIFitem(self.fp, '\n# PHASE INFO FOR HISTOGRAM '+hist)
4304                    # loop over phases, add a block header if there is more than one phase
4305                    for j,phasenam in enumerate(sorted(self.Phases.keys())):
4306                        pId = self.Phases[phasenam]['pId']
4307                        if hist not in self.Phases[phasenam]['Histograms']: continue
4308                        if not self.Phases[phasenam]['Histograms'][hist]['Use']: continue
4309                        WriteCIFitem(self.fp, '\n# phase info for '+str(phasenam) + ' follows')
4310                        if not phaseWithHist:
4311                            WriteCIFitem(self.fp, 'data_'+self.CIFname+"_hist"+str(i)+"_phase"+str(j))
4312                            WriteCIFitem(self.fp, '_pd_block_id',phaseBlockName[pId])
4313                            WriteCIFitem(self.fp, '')
4314
4315                        WriteSeqPhaseVals(phasenam,self.Phases[phasenam],pId,hist)
4316
4317                        # preferred orientation & profile terms
4318                        if self.ifPWDR:
4319                            #SH = FormatSH(phasenam)     # TODO: needs to use seqData
4320                            #MD = FormatHAPpo(phasenam)  # TODO: switch to seqData
4321                            #if SH and MD:
4322                            #    WriteCIFitem(self.fp, '_pd_proc_ls_pref_orient_corr', SH + '\n' + MD)
4323                            #elif SH or MD:
4324                            #    WriteCIFitem(self.fp, '_pd_proc_ls_pref_orient_corr', SH + MD)
4325                            #else:
4326                            #    WriteCIFitem(self.fp, '_pd_proc_ls_pref_orient_corr', 'none')
4327                            # report sample profile terms for all histograms with current phase
4328                            if phaseWithHist:
4329                                PP = FormatInstProfile(histblk["Instrument Parameters"],histblk['hId'])
4330                                PP += '\n'
4331                            else:
4332                                PP = ''
4333                            PP += FormatPhaseProfile(phasenam,hist)
4334                            WriteCIFitem(self.fp, '\n_pd_proc_ls_profile_function',PP)
4335            finally:
4336                dlg.Destroy()
4337        else:
4338            #======================================================================
4339            #### multiblock: multiple phases and/or histograms export
4340            #======================================================================
4341            oneblock = False
4342            for phasenam in sorted(self.Phases.keys()):
4343                rId = phasedict['ranId']
4344                if rId in self.CellHistSelection: continue
4345                self.CellHistSelection[rId] = self._CellSelectHist(phasenam)
4346            nsteps = 1 + len(self.Phases) + len(self.powderDict) + len(self.xtalDict)
4347            try:
4348                dlg = wx.ProgressDialog('CIF progress','starting',nsteps,parent=self.G2frame)
4349                dlg.CenterOnParent()
4350
4351                # publication info
4352                step = 1
4353                dlg.Update(step,"Exporting overall section")
4354                WriteCIFitem(self.fp, '\ndata_'+self.CIFname+'_publ')
4355                WriteAudit()
4356                WriteCIFitem(self.fp, '_pd_block_id',
4357                             str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
4358                             str(self.shortauthorname) + "|PubInfo")
4359                writeCIFtemplate(self.OverallParms['Controls'],'publ') #insert the publication template
4360                # ``template_publ.cif`` or a modified version
4361               
4362                # overall info -- it is not strictly necessary to separate this from the previous
4363                # publication block, but I think this makes sense
4364               
4365                WriteCIFitem(self.fp, 'data_'+str(self.CIFname)+'_overall')
4366                WriteCIFitem(self.fp, '_pd_block_id',
4367                             str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
4368                             str(self.shortauthorname) + "|Overall")
4369                if MM:
4370                    WriteOverallMM()
4371                else:
4372                    WriteOverall()
4373                #============================================================
4374                WriteCIFitem(self.fp, '# POINTERS TO PHASE AND/OR HISTOGRAM BLOCKS')
4375                datablockidDict = {} # save block names here -- N.B. check for conflicts between phase & hist names (unlikely!)
4376                # loop over phase blocks
4377                if len(self.Phases) > 1:
4378                    loopprefix = ''
4379                    WriteCIFitem(self.fp, 'loop_   _pd_phase_block_id')
4380                    for phasenam in sorted(self.Phases.keys()):
4381                        i = self.Phases[phasenam]['pId']
4382                        datablockidDict[phasenam] = (str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
4383                                                         'phase_'+ str(i) + '|' + str(self.shortauthorname))
4384                        WriteCIFitem(self.fp, loopprefix,datablockidDict[phasenam])
4385                else:    # phase in overall block
4386                    for phasenam in sorted(self.Phases.keys()): break
4387                    datablockidDict[phasenam] = (str(self.CIFdate) + "|" + str(self.CIFname) + "|" +
4388                                                     str(self.shortauthorname) + "|Overall")
4389                # loop over data blocks
4390                if len(self.powderDict) + len(self.xtalDict) > 1:
4391                    loopprefix = ''
4392                    WriteCIFitem(self.fp, 'loop_   _pd_block_diffractogram_id')
4393