source: trunk/GSASIIconstrGUI.py @ 4455

Last change on this file since 4455 was 4455, checked in by toby, 3 years ago

fix and speed up the "Complete Molecule" command; w/o for mult calc bug

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 158.8 KB
Line 
1# -*- coding: utf-8 -*-
2#GSASIIconstrGUI - constraint GUI routines
3########### SVN repository information ###################
4# $Date: 2020-06-02 00:16:10 +0000 (Tue, 02 Jun 2020) $
5# $Author: toby $
6# $Revision: 4455 $
7# $URL: trunk/GSASIIconstrGUI.py $
8# $Id: GSASIIconstrGUI.py 4455 2020-06-02 00:16:10Z toby $
9########### SVN repository information ###################
10'''
11*GSASIIconstrGUI: Constraint GUI routines*
12------------------------------------------
13
14Used to define constraints and rigid bodies.
15
16'''
17from __future__ import division, print_function
18import platform
19import sys
20import copy
21import os.path
22import wx
23import wx.grid as wg
24import wx.lib.scrolledpanel as wxscroll
25import wx.lib.gridmovers as gridmovers
26import random as ran
27import numpy as np
28import numpy.ma as ma
29import numpy.linalg as nl
30import GSASIIpath
31GSASIIpath.SetVersionNumber("$Revision: 4455 $")
32import GSASIIElem as G2elem
33import GSASIIElemGUI as G2elemGUI
34import GSASIIstrIO as G2stIO
35import GSASIImapvars as G2mv
36import GSASIImath as G2mth
37import GSASIIlattice as G2lat
38import GSASIIdataGUI as G2gd
39import GSASIIctrlGUI as G2G
40import GSASIIfiles as G2fl
41import GSASIIplot as G2plt
42import GSASIIobj as G2obj
43import GSASIIspc as G2spc
44import GSASIIpy3 as G2py3
45import GSASIIphsGUI as G2phG
46import GSASIIscriptable as G2sc
47VERY_LIGHT_GREY = wx.Colour(235,235,235)
48WACV = wx.ALIGN_CENTER_VERTICAL
49
50class G2BoolEditor(wg.GridCellBoolEditor):
51    '''Substitute for wx.grid.GridCellBoolEditor except toggles
52    grid items immediately when opened, updates grid & table contents after every
53    item change
54    '''
55    def __init__(self):
56        self.saveVals = None
57        wx.grid.GridCellBoolEditor.__init__(self)
58
59
60    def Create(self, parent, id, evtHandler):
61        '''Create the editing control (wx.CheckBox) when cell is opened
62        for edit
63        '''
64        self._tc = wx.CheckBox(parent, -1, "")
65        self._tc.Bind(wx.EVT_CHECKBOX, self.onCheckSet)
66        self.SetControl(self._tc)
67        if evtHandler:
68            self._tc.PushEventHandler(evtHandler)
69
70    def onCheckSet(self, event):
71        '''Callback used when checkbox is toggled.
72        Makes change to table immediately (creating event)
73        '''
74        if self.saveVals:
75            self.ApplyEdit(*self.saveVals)
76
77       
78    def SetSize(self, rect):
79        '''Set position/size the edit control within the cell's rectangle.
80        '''
81#        self._tc.SetDimensions(rect.x, rect.y, rect.width+2, rect.height+2, # older
82        self._tc.SetSize(rect.x, rect.y, rect.width+2, rect.height+2,
83                               wx.SIZE_ALLOW_MINUS_ONE)
84
85    def BeginEdit(self, row, col, grid):
86        '''Prepares the edit control by loading the initial
87        value from the table (toggles it since you would not
88        click on it if you were not planning to change it),
89        buts saves the original, pre-change value.
90        Makes change to table immediately.
91        Saves the info needed to make updates in self.saveVals.
92        Sets the focus.
93        '''
94        self.startValue = int(grid.GetTable().GetValue(row, col))
95        self.saveVals = row, col, grid
96        # invert state and set in editor
97        if self.startValue:
98            grid.GetTable().SetValue(row, col, 0)
99            self._tc.SetValue(0)
100        else:
101            grid.GetTable().SetValue(row, col, 1)
102            self._tc.SetValue(1)
103        self._tc.SetFocus()
104        self.ApplyEdit(*self.saveVals)
105
106    def EndEdit(self, row, col, grid, oldVal=None):
107        '''End editing the cell.  This is supposed to
108        return None if the value has not changed, but I am not
109        sure that actually works.
110        '''
111        val = int(self._tc.GetValue())
112        if val != oldVal:   #self.startValue:?
113            return val
114        else:
115            return None
116
117    def ApplyEdit(self, row, col, grid):
118        '''Save the value into the table, and create event.
119        Called after EndEdit(), BeginEdit and onCheckSet.
120        '''
121        val = int(self._tc.GetValue())
122        grid.GetTable().SetValue(row, col, val) # update the table
123
124    def Reset(self):
125        '''Reset the value in the control back to its starting value.
126        '''
127        self._tc.SetValue(self.startValue)
128
129    def StartingClick(self):
130        '''This seems to be needed for BeginEdit to work properly'''
131        pass
132
133    def Destroy(self):
134        "final cleanup"
135        super(G2BoolEditor, self).Destroy()
136
137    def Clone(self):
138        'required'
139        return G2BoolEditor()
140   
141class DragableRBGrid(wg.Grid):
142    '''Simple grid implentation for display of rigid body positions.
143
144    :param parent: frame or panel where grid will be placed
145    :param dict rb: dict with atom labels, types and positions
146    :param function onChange: a callback used every time a value in
147      rb is changed.
148    '''
149    def __init__(self, parent, rb, onChange=None):
150        #wg.Grid.__init__(self, parent, wx.ID_ANY,size=(-1,200))
151        wg.Grid.__init__(self, parent, wx.ID_ANY)
152        self.SetTable(RBDataTable(rb,onChange), True)
153        # Enable Row moving
154        gridmovers.GridRowMover(self)
155        self.Bind(gridmovers.EVT_GRID_ROW_MOVE, self.OnRowMove, self)
156        self.SetColSize(0, 60)
157        self.SetColSize(1, 40)
158        self.SetColSize(2, 35)
159        for r in range(len(rb['RBlbls'])):
160            self.SetReadOnly(r,0,isReadOnly=True)
161            self.SetCellEditor(r, 1, G2BoolEditor())           
162            self.SetCellRenderer(r, 1, wg.GridCellBoolRenderer())
163            self.SetReadOnly(r,2,isReadOnly=True)
164            self.SetCellEditor(r,3, wg.GridCellFloatEditor())
165            self.SetCellEditor(r,4, wg.GridCellFloatEditor())
166            self.SetCellEditor(r,6, wg.GridCellFloatEditor())
167
168    def OnRowMove(self,evt):
169        'called when a row move needs to take place'
170        frm = evt.GetMoveRow()          # Row being moved
171        to = evt.GetBeforeRow()         # Before which row to insert
172        self.GetTable().MoveRow(frm,to)
173       
174    def completeEdits(self):
175        'complete any outstanding edits'
176        if self.IsCellEditControlEnabled(): # complete any grid edits in progress
177            #if GSASIIpath.GetConfigValue('debug'): print ('Completing grid edit')
178            self.SaveEditControlValue()
179            self.HideCellEditControl()
180            self.DisableCellEditControl()
181           
182class RBDataTable(wg.GridTableBase):
183    '''A Table to support :class:`DragableRBGrid`
184    '''
185    def __init__(self,rb,onChange):
186        wg.GridTableBase.__init__(self)
187        self.colLabels = ['Label','Select','Type','x','y','z']
188        self.coords = rb['RBcoords']
189        self.labels = rb['RBlbls']
190        self.types = rb['RBtypes']
191        self.index = rb['RBindex']
192        self.select = rb['RBselection']
193        self.onChange = onChange
194
195    # required methods
196    def GetNumberRows(self):
197        return len(self.labels)
198    def GetNumberCols(self):
199        return len(self.colLabels)
200    def IsEmptyCell(self, row, col):
201        return False
202    def GetValue(self, row, col):
203        row = self.index[row]
204        if col == 0:
205            return self.labels[row]
206        elif col == 1:
207            return int(self.select[row])
208        elif col == 2:
209            return self.types[row]
210        else:
211            return '{:.5f}'.format(self.coords[row][col-3])
212    def SetValue(self, row, col, value):
213        row = self.index[row]
214        try:
215            if col == 0:
216                self.labels[row] = value
217            elif col == 1:
218                self.select[row] = int(value)
219            elif col == 2:
220                self.types[row] = value
221            else:
222                self.coords[row][col-3] = float(value)
223        except:
224            pass
225        if self.onChange:
226            self.onChange()
227    # Display column & row labels
228    def GetColLabelValue(self, col):
229        return self.colLabels[col]
230    def GetRowLabelValue(self,row):
231        return str(row+1)
232
233    # Implement "row movement" by updating the pointer array
234    def MoveRow(self,frm,to):
235        grid = self.GetView()
236        if grid:
237            move = self.index[frm]
238            del self.index[frm]
239            if frm > to:
240                self.index.insert(to,move)
241            else:
242                self.index.insert(to-1,move)
243           
244            # Notify the grid
245            grid.BeginBatch()
246            msg = wg.GridTableMessage(
247                    self, wg.GRIDTABLE_NOTIFY_ROWS_DELETED, frm, 1
248                    )
249            grid.ProcessTableMessage(msg)
250            msg = wg.GridTableMessage(
251                    self, wg.GRIDTABLE_NOTIFY_ROWS_INSERTED, to, 1
252                    )
253            grid.ProcessTableMessage(msg)
254            grid.EndBatch()
255        if self.onChange:
256            self.onChange()
257
258def MakeDrawAtom(data,atom):
259    'Convert atom to format needed to draw it'
260    generalData = data['General']
261    if generalData['Type'] in ['nuclear','faulted',]:
262        atomInfo = [atom[:2]+atom[3:6]+['1',]+['vdW balls',]+
263                    ['',]+[[255,255,255],]+atom[9:]+[[],[]]][0]
264    ct,cs = [1,8]         #type & color
265    atNum = generalData['AtomTypes'].index(atom[ct])
266    atomInfo[cs] = list(generalData['Color'][atNum])
267    return atomInfo
268
269class ConstraintDialog(wx.Dialog):
270    '''Window to edit Constraint values
271    '''
272    def __init__(self,parent,title,text,data,separator='*',varname="",varyflag=False):
273        wx.Dialog.__init__(self,parent,-1,'Edit '+title, 
274            pos=wx.DefaultPosition,style=wx.DEFAULT_DIALOG_STYLE)
275        self.data = data[:]
276        self.newvar = [varname,varyflag]
277        panel = wx.Panel(self)
278        mainSizer = wx.BoxSizer(wx.VERTICAL)
279        topLabl = wx.StaticText(panel,-1,text)
280        mainSizer.Add((10,10),1)
281        mainSizer.Add(topLabl,0,wx.ALIGN_CENTER_VERTICAL|wx.LEFT,10)
282        mainSizer.Add((10,10),1)
283        dataGridSizer = wx.FlexGridSizer(cols=3,hgap=2,vgap=2)
284        self.OkBtn = wx.Button(panel,wx.ID_OK)
285        for id in range(len(self.data)):
286            lbl1 = lbl = str(self.data[id][1])
287            if lbl[-1] != '=': lbl1 = lbl + ' ' + separator + ' '
288            name = wx.StaticText(panel,wx.ID_ANY,lbl1,
289                                 style=wx.ALIGN_RIGHT)
290            scale = G2G.ValidatedTxtCtrl(panel,self.data[id],0,
291                                          typeHint=float,
292                                          OKcontrol=self.DisableOK)
293            dataGridSizer.Add(name,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,5)
294            dataGridSizer.Add(scale,0,wx.RIGHT,3)
295            if ':' in lbl:
296                dataGridSizer.Add(
297                    wx.StaticText(panel,-1,G2obj.fmtVarDescr(lbl)),
298                    0,wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,3)
299            else:
300                dataGridSizer.Add((-1,-1))
301        if title == 'New Variable':
302            name = wx.StaticText(panel,wx.ID_ANY,"New variable's\nname (optional)",
303                                 style=wx.ALIGN_CENTER)
304            scale = G2G.ValidatedTxtCtrl(panel,self.newvar,0,
305                                          typeHint=str,notBlank=False)
306            dataGridSizer.Add(name,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,5)
307            dataGridSizer.Add(scale,0,wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,3)
308            self.refine = wx.CheckBox(panel,label='Refine?')
309            self.refine.SetValue(self.newvar[1]==True)
310            self.refine.Bind(wx.EVT_CHECKBOX, self.OnCheckBox)
311            dataGridSizer.Add(self.refine,0,wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,3)
312           
313        mainSizer.Add(dataGridSizer,0,wx.EXPAND)
314        self.OkBtn.Bind(wx.EVT_BUTTON, self.OnOk)
315        self.OkBtn.SetDefault()
316        cancelBtn = wx.Button(panel,wx.ID_CANCEL)
317        cancelBtn.Bind(wx.EVT_BUTTON, self.OnCancel)
318        btnSizer = wx.BoxSizer(wx.HORIZONTAL)
319        btnSizer.Add((20,20),1)
320        btnSizer.Add(self.OkBtn)
321        btnSizer.Add((20,20),1)
322        btnSizer.Add(cancelBtn)
323        btnSizer.Add((20,20),1)
324
325        mainSizer.Add(btnSizer,0,wx.EXPAND|wx.BOTTOM|wx.TOP, 10)
326        panel.SetSizer(mainSizer)
327        panel.Fit()
328        self.Fit()
329        self.CenterOnParent()
330       
331    def DisableOK(self,setting):
332        if setting:
333            self.OkBtn.Enable()
334        else:
335            self.OkBtn.Disable()
336
337    def OnCheckBox(self,event):
338        self.newvar[1] = self.refine.GetValue()
339
340    def OnOk(self,event):
341        parent = self.GetParent()
342        parent.Raise()
343        self.EndModal(wx.ID_OK)             
344
345    def OnCancel(self,event):
346        parent = self.GetParent()
347        parent.Raise()
348        self.EndModal(wx.ID_CANCEL)             
349
350    def GetData(self):
351        return self.data
352       
353################################################################################
354#####  Constraints
355################################################################################           
356def UpdateConstraints(G2frame,data):
357    '''Called when Constraints tree item is selected.
358    Displays the constraints in the data window
359    '''
360       
361    def FindEquivVarb(name,nameList):
362        'Creates a list of variables appropriate to constrain with name'
363        outList = []
364        #phlist = []
365        items = name.split(':')
366        namelist = [items[2],]
367        if 'dA' in name:
368            namelist = ['dAx','dAy','dAz']
369        elif 'AU' in name:
370            namelist = ['AUiso','AU11','AU22','AU33','AU12','AU13','AU23']
371        elif 'AM' in name:
372            namelist = ['AMx','AMy','AMz']
373        elif items[-1] in ['A0','A1','A2','A3','A4','A5']:
374            namelist = ['A0','A1','A2','A3','A4','A5']
375        elif items[-1] in ['D11','D22','D33','D12','D13','D23']:
376            namelist = ['D11','D22','D33','D12','D13','D23']
377        elif 'Tm' in name:
378            namelist = ['Tmin','Tmax']
379        elif 'MX' in name or 'MY' in name or 'MZ' in name:
380            namelist = ['MXcos','MYcos','MZcos','MXsin','MYsin','MZsin']
381        elif 'mV' in name:
382            namelist = ['mV0','mV1','mV2']
383        elif 'RB' in name:
384            rbfx = 'RB'+items[2][2]
385            if 'T' in name and 'Tr' not in name:
386                namelist = [rbfx+'T11',rbfx+'T22',rbfx+'T33',rbfx+'T12',rbfx+'T13',rbfx+'T23']
387            if 'L' in name:
388                namelist = [rbfx+'L11',rbfx+'L22',rbfx+'L33',rbfx+'L12',rbfx+'L13',rbfx+'L23']
389            if 'S' in name:
390                namelist = [rbfx+'S12',rbfx+'S13',rbfx+'S21',rbfx+'S23',rbfx+'S31',rbfx+'S32',rbfx+'SAA',rbfx+'SBB']
391            if 'U' in name:
392                namelist = [rbfx+'U',]
393
394        for item in nameList:
395            keys = item.split(':')
396            #if keys[0] not in phlist:
397            #    phlist.append(keys[0])
398            if items[1] == '*' and keys[2] in namelist: # wildcard -- select only sequential options
399                keys[1] = '*'
400                mitem = ':'.join(keys)
401                if mitem == name: continue
402                if mitem not in outList: outList.append(mitem)
403            elif keys[2] in namelist and item != name:
404                outList.append(item)
405        return outList
406       
407    def SelectVarbs(page,FrstVarb,varList,legend,constType):
408        '''Select variables used in constraints after one variable has
409        been selected. This routine determines the appropriate variables to be
410        used based on the one that has been selected and asks for more to be added.
411
412        It then creates the constraint and adds it to the constraints list.
413       
414        Called from OnAddEquivalence, OnAddFunction & OnAddConstraint (all but
415        OnAddHold)
416
417        :param list page: defines calling page -- type of variables to be used
418        :parm GSASIIobj.G2VarObj FrstVarb: reference to first selected variable
419        :param list varList: list of other appropriate variables to select from
420        :param str legend: header for selection dialog
421        :param str constType: type of constraint to be generated
422        :returns: a constraint, as defined in
423          :ref:`GSASIIobj <Constraint_definitions_table>`
424        '''
425        choices = [[i]+list(G2obj.VarDescr(i)) for i in varList]
426        meaning = G2obj.getDescr(FrstVarb.name)
427        if not meaning:
428            meaning = "(no definition found!)"
429        l = str(FrstVarb).split(':')
430        # make lists of phases & histograms to iterate over
431        phaselist = [l[0]]
432        if l[0]:
433            phaselbl = ['phase #'+l[0]]
434            if len(Phases) > 1:
435                phaselist += ['all'] 
436                phaselbl += ['all phases']
437        else:
438            phaselbl = ['']
439        histlist = [l[1]]
440        if l[1] == '*':
441            pass
442        elif l[1]:
443            histlbl = ['histogram #'+l[1]]
444            if len(Histograms) > 1:
445                histlist += ['all']
446                histlbl += ['all histograms']
447                typ = Histograms[G2obj.LookupHistName(l[1])[0]]['Instrument Parameters'][0]['Type'][1]
448                i = 0
449                for hist in Histograms:
450                    if Histograms[hist]['Instrument Parameters'][0]['Type'][1] == typ: i += 1
451                if i > 1:
452                    histlist += ['all='+typ]
453                    histlbl += ['all '+typ+' histograms']
454        else:
455            histlbl = ['']
456        # make a list of equivalent parameter names
457        nameList = [FrstVarb.name]
458        for var in varList:
459            nam = var.split(":")[2]
460            if nam not in nameList: nameList += [nam]
461        # add "wild-card" names to the list of variables
462        if l[1] == '*':
463            pass
464        elif page[1] == 'phs':
465            if 'RB' in FrstVarb.name:
466                pass
467            elif FrstVarb.atom is None:
468                for nam in nameList:
469                    for ph,plbl in zip(phaselist,phaselbl):
470                        if plbl: plbl = 'For ' + plbl
471                        var = ph+"::"+nam
472                        if var == str(FrstVarb) or var in varList: continue
473                        varList += [var]
474                        choices.append([var,plbl,meaning])
475            else:
476                for nam in nameList:
477                    for ph,plbl in zip(phaselist,phaselbl):
478                        if plbl: plbl = ' in ' + plbl
479                        for atype in ['']+TypeList:
480                            if atype:
481                                albl = "For "+atype+" atoms"
482                                akey = "all="+atype                       
483                            else:
484                                albl = "For all atoms"
485                                akey = "all"
486                            var = ph+"::"+nam+":"+akey
487                            if var == str(FrstVarb) or var in varList: continue
488                            varList += [var]
489                            choices.append([var,albl+plbl,meaning])
490        elif page[1] == 'hap':
491            if FrstVarb.name == "Scale":
492                meaning = "Phase fraction"
493            for nam in nameList:
494                for ph,plbl in zip(phaselist,phaselbl):
495                    if plbl: plbl = 'For ' + plbl
496                    for hst,hlbl in zip(histlist,histlbl):
497                        if hlbl:
498                            if plbl:
499                                hlbl = ' in ' + hlbl
500                            else:
501                                hlbl = 'For ' + hlbl                               
502                            var = ph+":"+hst+":"+nam
503                            if var == str(FrstVarb) or var in varList: continue
504                            varList += [var]
505                            choices.append([var,plbl+hlbl,meaning])
506        elif page[1] == 'hst':
507            if FrstVarb.name == "Scale":
508                meaning = "Scale factor"
509            for nam in nameList:
510                for hst,hlbl in zip(histlist,histlbl):
511                    if hlbl:
512                        hlbl = 'For ' + hlbl                               
513                        var = ":"+hst+":"+nam
514                        if var == str(FrstVarb) or var in varList: continue
515                        varList += [var]
516                        choices.append([var,hlbl,meaning])
517        elif page[1] == 'glb' or page[1] == 'sym':
518            pass
519        else:
520            raise Exception('Unknown constraint page '+ page[1])                   
521        if len(choices):
522            l1 = l2 = 1
523            for i1,i2,i3 in choices:
524                l1 = max(l1,len(i1))
525                l2 = max(l2,len(i2))
526            fmt = "{:"+str(l1)+"s} {:"+str(l2)+"s} {:s}"
527            atchoices = [fmt.format(*i1) for i1 in choices] # reformat list as str with columns
528            dlg = G2G.G2MultiChoiceDialog(
529                G2frame,legend,
530                'Constrain '+str(FrstVarb)+' with...',atchoices,
531                toggle=False,size=(625,400),monoFont=True)
532            dlg.CenterOnParent()
533            res = dlg.ShowModal()
534            Selections = dlg.GetSelections()[:]
535            dlg.Destroy()
536            if res != wx.ID_OK: return []
537            if len(Selections) == 0:
538                dlg = wx.MessageDialog(
539                    G2frame,
540                    'No variables were selected to include with '+str(FrstVarb),
541                    'No variables')
542                dlg.CenterOnParent()
543                dlg.ShowModal()
544                dlg.Destroy()
545                return []
546        else:
547            dlg = wx.MessageDialog(
548                G2frame,
549                'There are no appropriate variables to include with '+str(FrstVarb),
550                'No variables')
551            dlg.CenterOnParent()
552            dlg.ShowModal()
553            dlg.Destroy()
554            return []
555        # now process the variables provided by the user
556        varbs = [str(FrstVarb),] # list of selected variables
557        for sel in Selections:
558            var = varList[sel]
559            # phase(s) included
560            l = var.split(':')
561            if l[0] == "all":
562                phlist = [str(Phases[phase]['pId']) for phase in Phases]
563            else:
564                phlist = [l[0]]
565            # histogram(s) included
566            if l[1] == "all":
567                hstlist = [str(Histograms[hist]['hId']) for hist in Histograms]
568            elif '=' in l[1]:
569                htyp = l[1].split('=')[1]
570                hstlist = [str(Histograms[hist]['hId']) for hist in Histograms if 
571                           Histograms[hist]['Instrument Parameters'][0]['Type'][1] == htyp]
572            else:
573                hstlist = [l[1]]
574            if len(l) == 3:
575                for ph in phlist:
576                    for hst in hstlist:
577                        var = ph + ":" + hst + ":" + l[2]
578                        if var in varbs: continue
579                        varbs.append(var)
580            else: # constraints with atoms or rigid bodies
581                if len(l) == 5: # rigid body parameter
582                    var = ':'.join(l)
583                    if var in varbs: continue
584                    varbs.append(var)
585                elif l[3] == "all":
586                    for ph in phlist:
587                        key = G2obj.LookupPhaseName(ph)[0]
588                        for hst in hstlist: # should be blank
589                            for iatm,at in enumerate(Phases[key]['Atoms']):
590                                var = ph + ":" + hst + ":" + l[2] + ":" + str(iatm)
591                                if var in varbs: continue
592                                varbs.append(var)
593                elif '=' in l[3]:
594                    for ph in phlist:
595                        key = G2obj.LookupPhaseName(ph)[0]
596                        cx,ct,cs,cia = Phases[key]['General']['AtomPtrs']
597                        for hst in hstlist: # should be blank
598                            atyp = l[3].split('=')[1]
599                            for iatm,at in enumerate(Phases[key]['Atoms']):
600                                if at[ct] != atyp: continue
601                                var = ph + ":" + hst + ":" + l[2] + ":" + str(iatm)
602                                if var in varbs: continue
603                                varbs.append(var)
604                else:
605                    for ph in phlist:
606                        key = G2obj.LookupPhaseName(ph)[0]
607                        for hst in hstlist: # should be blank
608                            var = ph + ":" + hst + ":" + l[2] + ":" + l[3]
609                            if var in varbs: continue
610                            varbs.append(var)
611        if len(varbs) >= 1 or 'constraint' in constType:
612            constr = [[1.0,FrstVarb]]
613            for item in varbs[1:]:
614                constr += [[1.0,G2obj.G2VarObj(item)]]
615            if 'equivalence' in constType:
616                return [constr+[None,None,'e']]
617            elif 'function' in constType:
618                return [constr+[None,False,'f']]
619            elif 'constraint' in constType:
620                return [constr+[1.0,None,'c']]
621            else:
622                raise Exception('Unknown constraint type: '+str(constType))
623        else:
624            dlg = wx.MessageDialog(
625                G2frame,
626                'There are no selected variables to include with '+str(FrstVarb),
627                'No variables')
628            dlg.CenterOnParent()
629            dlg.ShowModal()
630            dlg.Destroy()
631        return []
632   
633    def ConstraintsLoad(data,newcons=[]):
634        '''Load all constraints. Constraints based on symmetry (etc.)
635        are generated by running :func:`GSASIIstrIO.GetPhaseData`.
636        '''
637        G2mv.InitVars()
638        #Find all constraints
639        constraintSet = []
640        for key in data:
641            if key.startswith('_'): continue
642            constraintSet += data[key]
643        if newcons:
644            constraintSet = constraintSet + newcons
645        constDictList,fixedList,ignored = G2stIO.ProcessConstraints(constraintSet)
646        # generate symmetry constraints to check for conflicts
647        rigidbodyDict = G2frame.GPXtree.GetItemPyData(   
648            G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Rigid bodies'))
649        rbIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[]})
650        rbVary,rbDict = G2stIO.GetRigidBodyModels(rigidbodyDict,Print=False)
651        Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables,MFtables,maxSSwave = G2stIO.GetPhaseData(
652            Phases,RestraintDict=None,rbIds=rbIds,Print=False) # generates atom symmetry constraints
653        return constDictList,phaseDict,fixedList
654           
655    def ConstraintsCheck(data,newcons=[]):
656        '''Load constraints & check them for errors. Since error checking
657        can cause changes in constraints in case of repairable conflicts
658        between equivalences, reload the constraints again after the check.
659        This could probably be done more effectively, but only reloading when
660        needed, but the reload should be quick.
661        '''
662        constDictList,phaseDict,fixedList = ConstraintsLoad(data,newcons)
663        msg = G2mv.EvaluateMultipliers(constDictList,phaseDict)
664        if msg:
665            return 'Unable to interpret multiplier(s): '+msg,''
666        res = G2mv.CheckConstraints('',constDictList,fixedList)
667        # reload constraints in case any were merged in MoveConfEquiv
668        ConstraintsLoad(data,newcons)
669        return res
670
671    def CheckAddedConstraint(newcons):
672        '''Check a new constraint that has just been input.
673        If there is an error display a message and discard the last entry
674
675        Since the varylist is not available, no warning messages
676        should be generated here
677
678        :returns: True if constraint should be added
679        '''
680       
681        errmsg,warnmsg = ConstraintsCheck(data,newcons)
682        if errmsg:
683            G2frame.ErrorDialog('Constraint Error',
684                'Error with newly added constraint:\n'+errmsg+
685                '\nIgnoring newly added constraint',parent=G2frame)
686            # reset error status
687            errmsg,warnmsg = ConstraintsCheck(data)
688            if errmsg:
689                print (errmsg)
690                print (G2mv.VarRemapShow([],True))
691            return False
692        elif warnmsg:
693            print ('Unexpected contraint warning:\n'+warnmsg)
694        return True
695
696    def CheckChangedConstraint():
697        '''Check all constraints after an edit has been made.
698        If there is an error display a message and reject the change.
699
700        Since the varylist is not available, no warning messages
701        should be generated.
702       
703        :returns: True if the edit should be retained
704        '''
705        errmsg,warnmsg = ConstraintsCheck(data)
706        if errmsg:
707            G2frame.ErrorDialog('Constraint Error',
708                'Error after editing constraint:\n'+errmsg+
709                '\nDiscarding last constraint edit',parent=G2frame)
710            # reset error status
711            errmsg,warnmsg = ConstraintsCheck(data)
712            if errmsg:
713                print (errmsg)
714                print (G2mv.VarRemapShow([],True))
715            return False
716        elif warnmsg:
717            print ('Unexpected contraint warning:\n'+warnmsg)
718        return True
719             
720    def PageSelection(page):
721        'Decode page reference'
722        if page[1] == "phs":
723            vartype = "phase"
724            varList = G2obj.removeNonRefined(phaseList)  # remove any non-refinable prms from list
725            constrDictEnt = 'Phase'
726        elif page[1] == "hap":
727            vartype = "Histogram*Phase"
728            varList = G2obj.removeNonRefined(hapList)  # remove any non-refinable prms from list
729            constrDictEnt = 'HAP'
730        elif page[1] == "hst":
731            vartype = "Histogram"
732            varList = G2obj.removeNonRefined(histList)  # remove any non-refinable prms from list
733            constrDictEnt = 'Hist'
734        elif page[1] == "glb":
735            vartype = "Global"
736            varList = G2obj.removeNonRefined(globalList)   # remove any non-refinable prms from list
737
738            constrDictEnt = 'Global'
739        elif page[1] == "sym":
740            return None,None,None
741        else:
742            raise Exception('Should not happen!')
743        return vartype,varList,constrDictEnt
744
745    def OnAddHold(event):
746        '''Create a new Hold constraint
747
748        Hold constraints allows the user to select one variable (the list of available
749        variables depends on which tab is currently active).
750        '''
751        page = G2frame.Page
752        vartype,varList,constrDictEnt = PageSelection(page)
753        if vartype is None: return
754        varList = G2obj.SortVariables(varList)
755        title1 = "Hold "+vartype+" variable"
756        if not varList:
757            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
758                parent=G2frame)
759            return
760        l2 = l1 = 1
761        for i in varList:
762            l1 = max(l1,len(i))
763            loc,desc = G2obj.VarDescr(i)
764            l2 = max(l2,len(loc))
765        fmt = "{:"+str(l1)+"s} {:"+str(l2)+"s} {:s}"
766        varListlbl = [fmt.format(i,*G2obj.VarDescr(i)) for i in varList]
767        #varListlbl = ["("+i+") "+G2obj.fmtVarDescr(i) for i in varList]
768        legend = "Select variables to hold (Will not be varied, even if vary flag is set)"
769        dlg = G2G.G2MultiChoiceDialog(G2frame,
770            legend,title1,varListlbl,toggle=False,size=(625,400),monoFont=True)
771        dlg.CenterOnParent()
772        if dlg.ShowModal() == wx.ID_OK:
773            for sel in dlg.GetSelections():
774                Varb = varList[sel]
775                VarObj = G2obj.G2VarObj(Varb)
776                newcons = [[[0.0,VarObj],None,None,'h']]
777                if CheckAddedConstraint(newcons):
778                    data[constrDictEnt] += newcons
779        dlg.Destroy()
780        wx.CallAfter(OnPageChanged,None)
781       
782    def OnAddEquivalence(event):
783        '''add an Equivalence constraint'''
784        page = G2frame.Page
785        vartype,varList,constrDictEnt = PageSelection(page)
786        if vartype is None: return
787        title1 = "Create equivalence constraint between "+vartype+" variables"
788        title2 = "Select additional "+vartype+" variable(s) to be equivalent with "
789        if not varList:
790            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
791                parent=G2frame)
792            return
793#        legend = "Select variables to make equivalent (only one of the variables will be varied when all are set to be varied)"
794        GetAddVars(page,title1,title2,varList,constrDictEnt,'equivalence')
795       
796    def OnAddAtomEquiv(event):
797        ''' Add equivalences between all parameters on atoms '''
798        page = G2frame.Page
799        vartype,varList,constrDictEnt = PageSelection(page)
800        if vartype is None: return
801        title1 = "Setup equivalent atom variables"
802        title2 = "Select additional atoms(s) to be equivalent with "
803        if not varList:
804            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
805                parent=G2frame)
806            return
807#        legend = "Select atoms to make equivalent (only one of the atom variables will be varied when all are set to be varied)"
808        GetAddAtomVars(page,title1,title2,varList,constrDictEnt,'equivalence')
809       
810    def OnAddRiding(event):
811        ''' Add riding equivalences between all parameters on atoms '''
812        page = G2frame.Page
813        vartype,varList,constrDictEnt = PageSelection(page)
814        if vartype is None: return
815        title1 = "Setup riding atoms "
816        title2 = "Select additional atoms(s) to ride on "
817        if not varList:
818            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
819                parent=G2frame)
820            return
821#        legend = "Select atoms to ride (only one of the atom variables will be varied when all are set to be varied)"
822        GetAddAtomVars(page,title1,title2,varList,constrDictEnt,'riding')
823   
824    def OnAddFunction(event):
825        '''add a Function (new variable) constraint'''
826        page = G2frame.Page
827        vartype,varList,constrDictEnt = PageSelection(page)
828        if vartype is None: return
829        title1 = "Setup new variable based on "+vartype+" variables"
830        title2 = "Include additional "+vartype+" variable(s) to be included with "
831        if not varList:
832            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
833                parent=G2frame)
834            return
835#        legend = "Select variables to include in a new variable (the new variable will be varied when all included variables are varied)"
836        GetAddVars(page,title1,title2,varList,constrDictEnt,'function')
837                       
838    def OnAddConstraint(event):
839        '''add a constraint equation to the constraints list'''
840        page = G2frame.Page
841        vartype,varList,constrDictEnt = PageSelection(page)
842        if vartype is None: return
843        title1 = "Creating constraint on "+vartype+" variables"
844        title2 = "Select additional "+vartype+" variable(s) to include in constraint with "
845        if not varList:
846            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
847                parent=G2frame)
848            return
849#        legend = "Select variables to include in a constraint equation (the values will be constrainted to equal a specified constant)"
850        GetAddVars(page,title1,title2,varList,constrDictEnt,'constraint')
851
852    def GetAddVars(page,title1,title2,varList,constrDictEnt,constType):
853        '''Get the variables to be added for OnAddEquivalence, OnAddFunction,
854        and OnAddConstraint. Then create and check the constraint.
855        '''
856        #varListlbl = ["("+i+") "+G2obj.fmtVarDescr(i) for i in varList]
857        if constType == 'equivalence':
858            omitVars = G2mv.GetDependentVars()
859        else:
860            omitVars = []
861        varList = G2obj.SortVariables([i for i in varList if i not in omitVars])
862        l2 = l1 = 1
863        for i in varList:
864            l1 = max(l1,len(i))
865            loc,desc = G2obj.VarDescr(i)
866            l2 = max(l2,len(loc))
867        fmt = "{:"+str(l1)+"s} {:"+str(l2)+"s} {:s}"
868        varListlbl = [fmt.format(i,*G2obj.VarDescr(i)) for i in varList]
869        dlg = G2G.G2SingleChoiceDialog(G2frame,'Select 1st variable:',
870            title1,varListlbl,monoFont=True,size=(625,400))
871        dlg.CenterOnParent()
872        if dlg.ShowModal() == wx.ID_OK:
873            if constType == 'equivalence':
874                omitVars = G2mv.GetDependentVars() + G2mv.GetIndependentVars()
875            sel = dlg.GetSelection()
876            FrstVarb = varList[sel]
877            VarObj = G2obj.G2VarObj(FrstVarb)
878            moreVarb = G2obj.SortVariables(FindEquivVarb(FrstVarb,[i for i in varList if i not in omitVars]))
879            newcons = SelectVarbs(page,VarObj,moreVarb,title2+FrstVarb,constType)
880            if len(newcons) > 0:
881                if CheckAddedConstraint(newcons):
882                    data[constrDictEnt] += newcons
883        dlg.Destroy()
884        wx.CallAfter(OnPageChanged,None)
885                       
886    def FindNeighbors(phase,FrstName,AtNames):
887        General = phase['General']
888        cx,ct,cs,cia = General['AtomPtrs']
889        Atoms = phase['Atoms']
890        atNames = [atom[ct-1] for atom in Atoms]
891        Cell = General['Cell'][1:7]
892        Amat,Bmat = G2lat.cell2AB(Cell)
893        atTypes = General['AtomTypes']
894        Radii = np.array(General['BondRadii'])
895        AtInfo = dict(zip(atTypes,Radii)) #or General['BondRadii']
896        Orig = atNames.index(FrstName.split()[1])
897        OType = Atoms[Orig][ct]
898        XYZ = G2mth.getAtomXYZ(Atoms,cx)       
899        Neigh = []
900        Dx = np.inner(Amat,XYZ-XYZ[Orig]).T
901        dist = np.sqrt(np.sum(Dx**2,axis=1))
902        sumR = AtInfo[OType]+0.5    #H-atoms only!
903        IndB = ma.nonzero(ma.masked_greater(dist-0.85*sumR,0.))
904        for j in IndB[0]:
905            if j != Orig:
906                Neigh.append(AtNames[j])
907        return Neigh
908       
909    def GetAddAtomVars(page,title1,title2,varList,constrDictEnt,constType):
910        '''Get the atom variables to be added for OnAddAtomEquiv. Then create and
911        check the constraints. Riding for H atoms only.
912        '''
913        Atoms = {G2obj.VarDescr(i)[0]:[] for i in varList if 'Atom' in G2obj.VarDescr(i)[0]}
914        for item in varList:
915            atName = G2obj.VarDescr(item)[0]
916            if atName in Atoms:
917                Atoms[atName].append(item)
918        AtNames = list(Atoms.keys())
919        AtNames.sort()
920        dlg = G2G.G2SingleChoiceDialog(G2frame,'Select 1st atom:',
921            title1,AtNames,monoFont=True,size=(625,400))
922        dlg.CenterOnParent()
923        FrstAtom = ''
924        if dlg.ShowModal() == wx.ID_OK:
925            sel = dlg.GetSelection()
926            FrstAtom = AtNames[sel]
927            if 'riding' in constType:
928                phaseName = (FrstAtom.split(' in ')[1]).strip()
929                phase = Phases[phaseName]
930                AtNames = FindNeighbors(phase,FrstAtom,AtNames)
931            else:
932                AtNames.remove(FrstAtom)
933        dlg.Destroy()
934        if FrstAtom == '':
935            print ('no atom selected')
936            return
937        dlg = G2G.G2MultiChoiceDialog(
938            G2frame,title2+FrstAtom,
939            'Constrain '+str(FrstAtom)+' with...',AtNames,
940            toggle=False,size=(625,400),monoFont=True)
941        if dlg.ShowModal() == wx.ID_OK:
942            Selections = dlg.GetSelections()[:]
943        else:
944            print ('no target atom selected')
945            dlg.Destroy()
946            return
947        dlg.Destroy()
948        for name in Atoms[FrstAtom]:
949            newcons = []
950            constr = []
951            if 'riding' in constType:
952                if 'AUiso' in name:
953                    constr = [[1.0,G2obj.G2VarObj(name)]]
954                elif 'AU11' in name:
955                    pass
956                elif 'AU' not in name:
957                    constr = [[1.0,G2obj.G2VarObj(name)]]
958            else:
959                constr = [[1.0,G2obj.G2VarObj(name)]]
960            pref = ':'+name.rsplit(':',1)[0].split(':',1)[1]    #get stuff between phase id & atom id
961            for sel in Selections:
962                name2 = Atoms[AtNames[sel]][0]
963                pid = name2.split(':',1)[0]                     #get phase id for 2nd atom
964                id = name2.rsplit(':',1)[-1]                    #get atom no. for 2nd atom
965                if 'riding' in constType:
966                    pref = pid+pref
967                    if 'AUiso' in pref:
968                        parts = pref.split('AUiso')
969                        constr += [[1.2,G2obj.G2VarObj('%s:%s'%(parts[0]+'AUiso',id))]]
970                    elif 'AU' not in pref:
971                        constr += [[1.0,G2obj.G2VarObj('%s:%s'%(pref,id))]]
972                else:
973                    constr += [[1.0,G2obj.G2VarObj('%s:%s'%(pid+pref,id))]]
974            if not constr:
975                continue
976            if 'frac' in pref and 'riding' not in constType:
977                newcons = [constr+[1.0,None,'c']]
978            else:
979                newcons = [constr+[None,None,'e']]
980            if len(newcons) > 0:
981                if CheckAddedConstraint(newcons):
982                    data[constrDictEnt] += newcons
983        wx.CallAfter(OnPageChanged,None)
984                       
985    def MakeConstraintsSizer(name,pageDisplay):
986        '''Creates a sizer displaying all of the constraints entered of
987        the specified type.
988
989        :param str name: the type of constraints to be displayed ('HAP',
990          'Hist', 'Phase', 'Global', 'Sym-Generated')
991        :param wx.Panel pageDisplay: parent panel for sizer
992        :returns: wx.Sizer created by method
993        '''
994        if name == 'Sym-Generated':         #show symmetry generated constraints
995            Sizer1 =  wx.BoxSizer(wx.VERTICAL)
996            Sizer1.Add(wx.StaticText(pageDisplay,wx.ID_ANY,
997                                    'Equivalences generated based on cell/space group input'))
998            Sizer1.Add((-1,5))
999            Sizer = wx.FlexGridSizer(0,2,0,0)
1000            Sizer1.Add(Sizer)
1001            for sym in G2mv.GetSymEquiv():
1002                Sizer.Add(wx.StaticText(pageDisplay,wx.ID_ANY,'EQUIV'),
1003                           0,wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER|wx.RIGHT|wx.LEFT,3)
1004                Sizer.Add(wx.StaticText(pageDisplay,wx.ID_ANY,sym))
1005                Sizer.Add((-1,-1))
1006                Sizer.Add((-1,2))
1007            return Sizer1
1008        constSizer = wx.FlexGridSizer(0,6,0,0)
1009        maxlen = 50 # characters before wrapping a constraint
1010        for Id,item in enumerate(data[name]):
1011            refineflag = False
1012            helptext = ""
1013            eqString = ['',]
1014            problemItem = False
1015            for term in item[:-3]:
1016                if str(term[1]) in G2mv.problemVars:
1017                    problemItem = True
1018            if item[-1] == 'h': # Hold on variable
1019                constSizer.Add((-1,-1),0)              # blank space for edit button
1020                typeString = 'FIXED'
1021                var = str(item[0][1])
1022                varMean = G2obj.fmtVarDescr(var)
1023                eqString[-1] =  var +'   '
1024                helptext = "Prevents variable:\n"+ var + " ("+ varMean + ")\nfrom being changed"
1025            elif item[-1] == 'f' or item[-1] == 'e' or item[-1] == 'c': # not true on original-style (2011?) constraints
1026                constEdit = wx.Button(pageDisplay,wx.ID_ANY,'Edit',style=wx.BU_EXACTFIT)
1027                constEdit.Bind(wx.EVT_BUTTON,OnConstEdit)
1028                constSizer.Add(constEdit,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)            # edit button
1029                Indx[constEdit.GetId()] = [Id,name]
1030                if item[-1] == 'f':
1031                    helptext = "A new variable"
1032                    if item[-3]:
1033                        helptext += " named "+str(item[-3])
1034                    helptext += " is created from a linear combination of the following variables:\n"
1035                    for term in item[:-3]:
1036                        var = str(term[1])
1037                        if len(eqString[-1]) > maxlen:
1038                            eqString.append(' ')
1039                        m = term[0]
1040                        if eqString[-1] != '':
1041                            if m >= 0:
1042                                eqString[-1] += ' + '
1043                            else:
1044                                eqString[-1] += ' - '
1045                                m = abs(m)
1046                        if m == 1:
1047                            eqString[-1] += '{:} '.format(var)
1048                        else:
1049                            eqString[-1] += '{:g}*{:} '.format(m,var)
1050                        varMean = G2obj.fmtVarDescr(var)
1051                        helptext += "\n" + var + " ("+ varMean + ")"
1052                    # Add ISODISTORT help items
1053                    if '_Explain' in data:
1054                        # this ignores the phase number. TODO: refactor that in
1055                        hlptxt = None
1056                        try:
1057                            hlptxt = data['_Explain'].get(item[-3])
1058                        except TypeError:
1059                            # note that phase RanId is item[-3].phase
1060                            hlptxt = data['_Explain'].get(str(item[-3].phase)+item[-3].name)
1061                        if hlptxt:
1062                            helptext += '\n\n'
1063                            helptext += hlptxt
1064                    # typeString = 'NEWVAR'
1065                    # if item[-3]:
1066                    #     eqString[-1] += ' = '+item[-3]
1067                    # else:
1068                    #     eqString[-1] += ' = New Variable'
1069                    if item[-3]:
1070                        typeString = str(item[-3]) + ' = '
1071                    else:
1072                        typeString = 'New Variable = '
1073                    #print 'refine',item[-2]
1074                    refineflag = True
1075                elif item[-1] == 'c':
1076                    helptext = "The following variables constrained to equal a constant:"
1077                    for term in item[:-3]:
1078                        var = str(term[1])
1079                        if len(eqString[-1]) > maxlen:
1080                            eqString.append(' ')
1081                        m = term[0]
1082                        if eqString[-1] != '':
1083                            if term[0] > 0:
1084                                eqString[-1] += ' + '
1085                            else:
1086                                eqString[-1] += ' - '
1087                                m = -term[0]
1088                        if m == 1:
1089                            eqString[-1] += '{:} '.format(var)
1090                        else:
1091                            eqString[-1] += '{:g}*{:} '.format(m,var)
1092                        varMean = G2obj.fmtVarDescr(var)
1093                        helptext += "\n" + var + " ("+ varMean + ")"
1094                    typeString = 'CONST'
1095                    eqString[-1] += ' = '+str(item[-3])
1096                elif item[-1] == 'e':
1097                    helptext = "The following variables are set to be equivalent, noting multipliers:"
1098                    normval = item[:-3][1][0]
1099                    for i,term in enumerate(item[:-3]):
1100                        var = str(term[1])
1101                        if term[0] == 0: term[0] = 1.0
1102                        if len(eqString[-1]) > maxlen:
1103                            eqString.append(' ')
1104                        if i == 0: # move independent variable to end
1105                            indepterm = term
1106                            continue
1107                        elif eqString[-1] != '':
1108                            eqString[-1] += ' = '
1109                        if normval/term[0] == 1:
1110                            eqString[-1] += '{:}'.format(var)
1111                        else:
1112                            eqString[-1] += '{:g}*{:} '.format(normval/term[0],var)
1113                        varMean = G2obj.fmtVarDescr(var)
1114                        helptext += "\n" + var + " ("+ varMean + ")"
1115                    if normval/indepterm[0] == 1:
1116                        eqString[-1] += ' = {:} '.format(indepterm[1])
1117                    else:
1118                        eqString[-1] += ' = {:g}*{:} '.format(normval/indepterm[0],str(indepterm[1]))
1119                    typeString = 'EQUIV'
1120                else:
1121                    print ('Unexpected constraint'+item)
1122               
1123            else:
1124                print ('Removing old-style constraints')
1125                data[name] = []
1126                return constSizer
1127            constDel = wx.Button(pageDisplay,wx.ID_ANY,'Delete',style=wx.BU_EXACTFIT)
1128            constDel.Bind(wx.EVT_BUTTON,OnConstDel)
1129            Indx[constDel.GetId()] = [Id,name]
1130            constSizer.Add(constDel,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)             # delete button
1131            if helptext:
1132                ch = G2G.HelpButton(pageDisplay,helptext)
1133                constSizer.Add(ch,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)
1134            else:
1135                constSizer.Add((-1,-1))
1136            if refineflag:
1137                ch = G2G.G2CheckBox(pageDisplay,'',item,-2)
1138                constSizer.Add(ch,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)
1139            else:
1140                constSizer.Add((-1,-1))               
1141            constSizer.Add(wx.StaticText(pageDisplay,wx.ID_ANY,typeString),
1142                           0,wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER|wx.RIGHT|wx.LEFT,3)
1143            if problemItem: eqString[-1] += ' -- Conflict: see console'
1144            if len(eqString) > 1:
1145                Eq = wx.BoxSizer(wx.VERTICAL)
1146                for s in eqString:
1147                    line = wx.StaticText(pageDisplay,wx.ID_ANY,s)
1148                    if problemItem: line.SetBackgroundColour(wx.YELLOW)
1149                    Eq.Add(line,0,wx.ALIGN_CENTER_VERTICAL)
1150                Eq.Add((-1,4))
1151            else:
1152                Eq = wx.StaticText(pageDisplay,wx.ID_ANY,eqString[0])
1153                if problemItem: Eq.SetBackgroundColour(wx.YELLOW)
1154            constSizer.Add(Eq,1,wx.ALIGN_CENTER_VERTICAL)
1155        return constSizer
1156               
1157    def OnConstDel(event):
1158        'Delete a constraint'
1159        Obj = event.GetEventObject()
1160        Id,name = Indx[Obj.GetId()]
1161        del data[name][Id]
1162        ConstraintsLoad(data)
1163        wx.CallAfter(OnPageChanged,None)
1164
1165    def OnConstEdit(event):
1166        '''Called to edit an individual contraint in response to a
1167        click on its Edit button
1168        '''
1169        Obj = event.GetEventObject()
1170        Id,name = Indx[Obj.GetId()]
1171        if data[name][Id][-1] == 'f':
1172            items = data[name][Id][:-3]
1173            constType = 'New Variable'
1174            if data[name][Id][-3]:
1175                varname = str(data[name][Id][-3])
1176            else:
1177                varname = ""
1178            lbl = 'Enter value for each term in constraint; sum = new variable'
1179            dlg = ConstraintDialog(G2frame,constType,lbl,items,
1180                                   varname=varname,varyflag=data[name][Id][-2])
1181        elif data[name][Id][-1] == 'c':
1182            items = data[name][Id][:-3]+[
1183                [str(data[name][Id][-3]),'fixed value =']]
1184            constType = 'Constraint'
1185            lbl = 'Edit value for each term in constant constraint sum'
1186            dlg = ConstraintDialog(G2frame,constType,lbl,items)
1187        elif data[name][Id][-1] == 'e':
1188            items = data[name][Id][:-3]
1189            constType = 'Equivalence'
1190            lbl = 'The following terms are set to be equal:'
1191            dlg = ConstraintDialog(G2frame,constType,lbl,items,'/')
1192        else:
1193            return
1194        try:
1195            prev = copy.deepcopy(data[name][Id])
1196            if dlg.ShowModal() == wx.ID_OK:
1197                result = dlg.GetData()
1198                for i in range(len(data[name][Id][:-3])):
1199                    if type(data[name][Id][i]) is tuple: # fix non-mutable construct
1200                        data[name][Id][i] = list(data[name][Id][i])
1201                    data[name][Id][i][0] = result[i][0]
1202                if data[name][Id][-1] == 'c':
1203                    data[name][Id][-3] = str(result[-1][0])
1204                elif data[name][Id][-1] == 'f':
1205                    data[name][Id][-2] = dlg.newvar[1]
1206                    if type(data[name][Id][-3]) is str:
1207                        # process the variable name to put in global form (::var)
1208                        varname = str(dlg.newvar[0]).strip().replace(' ','_')
1209                        if varname.startswith('::'):
1210                            varname = varname[2:]
1211                        varname = varname.replace(':',';')
1212                        if varname:
1213                            data[name][Id][-3] = varname
1214                        else:
1215                            data[name][Id][-3] = ''
1216                if not CheckChangedConstraint():
1217                    data[name][Id] = prev
1218            else:
1219                data[name][Id] = prev
1220        except:
1221            import traceback
1222            print (traceback.format_exc())
1223        finally:
1224            dlg.Destroy()
1225        wx.CallAfter(OnPageChanged,None)
1226   
1227    def UpdateConstraintPanel(panel,typ):
1228        '''Update the contents of the selected Constraint
1229        notebook tab. Called in :func:`OnPageChanged`
1230        '''
1231        if panel.GetSizer(): panel.GetSizer().Clear(True)
1232        Siz = wx.BoxSizer(wx.VERTICAL)
1233        Siz.Add((5,5),0)
1234        Siz.Add(MakeConstraintsSizer(typ,panel),1,wx.EXPAND)
1235        panel.SetSizer(Siz,True)
1236        Size = Siz.GetMinSize()
1237        Size[0] += 40
1238        Size[1] = max(Size[1],450) + 20
1239        panel.SetSize(Size)
1240        panel.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
1241        panel.Show()
1242
1243    def OnPageChanged(event):
1244        '''Called when a tab is pressed or when a "select tab" menu button is
1245        used (see RaisePage), or to refresh the current tab contents (event=None)
1246        '''
1247        if event:       #page change event!
1248            page = event.GetSelection()
1249        else: # called directly, get current page
1250            page = G2frame.constr.GetSelection()
1251        G2frame.constr.ChangeSelection(page)
1252        text = G2frame.constr.GetPageText(page)
1253        G2frame.dataWindow.ConstraintEdit.Enable(G2G.wxID_EQUIVALANCEATOMS,False)
1254#        G2frame.dataWindow.ConstraintEdit.Enable(G2G.wxID_ADDRIDING,False)
1255        if text == 'Histogram/Phase':
1256            enableEditCons = [False]+4*[True]
1257            G2frame.Page = [page,'hap']
1258            UpdateConstraintPanel(HAPConstr,'HAP')
1259        elif text == 'Histogram':
1260            enableEditCons = [False]+4*[True]
1261            G2frame.Page = [page,'hst']
1262            UpdateConstraintPanel(HistConstr,'Hist')
1263        elif text == 'Phase':
1264            enableEditCons = 5*[True]
1265            G2frame.Page = [page,'phs']
1266            G2frame.dataWindow.ConstraintEdit.Enable(G2G.wxID_EQUIVALANCEATOMS,True)
1267#            G2frame.dataWindow.ConstraintEdit.Enable(G2G.wxID_ADDRIDING,True)
1268            if 'DELETED' in str(PhaseConstr):   #seems to be no other way to do this (wx bug)
1269                if GSASIIpath.GetConfigValue('debug'):
1270                    print ('DBG_wx error: PhaseConstr not cleanly deleted after Refine')
1271                return
1272            UpdateConstraintPanel(PhaseConstr,'Phase')
1273        elif text == 'Global':
1274            enableEditCons = [False]+4*[True]
1275            G2frame.Page = [page,'glb']
1276            UpdateConstraintPanel(GlobalConstr,'Global')
1277        else:
1278            enableEditCons = 5*[False]
1279            G2frame.Page = [page,'sym']
1280            UpdateConstraintPanel(SymConstr,'Sym-Generated')
1281        # remove menu items when not allowed
1282        for obj,flag in zip(G2frame.dataWindow.ConstraintEdit.GetMenuItems(),enableEditCons): 
1283            obj.Enable(flag)
1284        G2frame.dataWindow.SetDataSize()
1285
1286    def RaisePage(event):
1287        'Respond to a "select tab" menu button'
1288        try:
1289            i = (G2G.wxID_CONSPHASE,
1290                 G2G.wxID_CONSHAP,
1291                 G2G.wxID_CONSHIST,
1292                 G2G.wxID_CONSGLOBAL,
1293                 G2G.wxID_CONSSYM,
1294                ).index(event.GetId())
1295            G2frame.constr.SetSelection(i)
1296            wx.CallAfter(OnPageChanged,None)
1297        except ValueError:
1298            print('Unexpected event in RaisePage')
1299
1300    def SetStatusLine(text):
1301        G2frame.GetStatusBar().SetStatusText(text,1)
1302
1303    # UpdateConstraints execution starts here
1304    if not data:
1305        data.update({'Hist':[],'HAP':[],'Phase':[],'Global':[]})       #empty dict - fill it
1306    if 'Global' not in data:                                            #patch
1307        data['Global'] = []
1308    # DEBUG code ########################################
1309    #import GSASIIconstrGUI
1310    #reload(GSASIIconstrGUI)
1311    #reload(G2obj)
1312    #reload(G2stIO)
1313    #import GSASIIstrMain
1314    #reload(GSASIIstrMain)   
1315    #reload(G2mv)
1316    #reload(G2gd)
1317    ###################################################
1318    Histograms,Phases = G2frame.GetUsedHistogramsAndPhasesfromTree()
1319    if not len(Phases) or not len(Histograms):
1320        dlg = wx.MessageDialog(G2frame,'You need both phases and histograms to see Constraints',
1321            'No phases or histograms')
1322        dlg.CenterOnParent()
1323        dlg.ShowModal()
1324        dlg.Destroy()
1325        return
1326    G2obj.IndexAllIds(Histograms,Phases)
1327    for p in Phases:
1328        if 'ISODISTORT' in Phases[p]:
1329            G2frame.dataWindow.ConstraintEdit.Enable(G2G.wxID_SHOWISO,True)
1330            break
1331    ##################################################################################
1332    # patch: convert old-style (str) variables in constraints to G2VarObj objects
1333    for key,value in data.items():
1334        if key.startswith('_'): continue
1335        j = 0
1336        for cons in value:
1337            #print cons             # DEBUG
1338            for i in range(len(cons[:-3])):
1339                if type(cons[i][1]) is str:
1340                    cons[i][1] = G2obj.G2VarObj(cons[i][1])
1341                    j += 1
1342        if j:
1343            print (str(key) + ': '+str(j)+' variable(s) as strings converted to objects')
1344    ##################################################################################
1345    rigidbodyDict = G2frame.GPXtree.GetItemPyData(
1346        G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Rigid bodies'))
1347    rbIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[]})
1348    rbVary,rbDict = G2stIO.GetRigidBodyModels(rigidbodyDict,Print=False)
1349    badPhaseParms = ['Ax','Ay','Az','Amul','AI/A','Atype','SHorder','AwaveType','FwaveType','PwaveType','MwaveType','Vol','isMag',]
1350    globalList = list(rbDict.keys())
1351    globalList.sort()
1352    try:
1353        AtomDict = dict([Phases[phase]['pId'],Phases[phase]['Atoms']] for phase in Phases)
1354    except KeyError:
1355        G2frame.ErrorDialog('Constraint Error','Constraints cannot be set until a cycle of least squares'+
1356                            ' has been run.\nWe suggest you refine a scale factor.')
1357        return
1358
1359    # create a list of the phase variables
1360    Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtable,BLtable,MFtable,maxSSwave = G2stIO.GetPhaseData(Phases,rbIds=rbIds,Print=False)
1361    phaseList = []
1362    for item in phaseDict:
1363        if item.split(':')[2] not in badPhaseParms:
1364            phaseList.append(item)
1365    phaseList.sort()
1366    phaseAtNames = {}
1367    phaseAtTypes = {}
1368    TypeList = []
1369    for item in phaseList:
1370        Split = item.split(':')
1371        if Split[2][:2] in ['AU','Af','dA','AM']:
1372            Id = int(Split[0])
1373            phaseAtNames[item] = AtomDict[Id][int(Split[3])][0]
1374            phaseAtTypes[item] = AtomDict[Id][int(Split[3])][1]
1375            if phaseAtTypes[item] not in TypeList:
1376                TypeList.append(phaseAtTypes[item])
1377        else:
1378            phaseAtNames[item] = ''
1379            phaseAtTypes[item] = ''
1380             
1381    # create a list of the hist*phase variables
1382    seqList = G2frame.testSeqRefineMode()
1383    if seqList: # for sequential refinement, only process 1st histgram in list
1384        histDict = {seqList[0]:Histograms[seqList[0]]}
1385    else:
1386        histDict = Histograms
1387    hapVary,hapDict,controlDict = G2stIO.GetHistogramPhaseData(Phases,histDict,Print=False,resetRefList=False)
1388    hapList = sorted([i for i in hapDict.keys() if i.split(':')[2] not in ('Type',)])
1389    if seqList: # convert histogram # to wildcard
1390        wildList = [] # list of variables with "*" for histogram number
1391        for i in hapList:
1392            s = i.split(':')
1393            if s[1] == "": continue
1394            s[1] = '*'
1395            sj = ':'.join(s)
1396            if sj not in wildList: wildList.append(sj)
1397        hapList = wildList
1398    histVary,histDict,controlDict = G2stIO.GetHistogramData(histDict,Print=False)
1399    histList = list(histDict.keys())
1400    histList.sort()
1401    if seqList: # convert histogram # to wildcard
1402        wildList = [] # list of variables with "*" for histogram number
1403        for i in histList:
1404            s = i.split(':')
1405            if s[1] == "": continue
1406            s[1] = '*'
1407            sj = ':'.join(s)
1408            if sj not in wildList: wildList.append(sj)
1409        histList = wildList       
1410    Indx = {}
1411    G2frame.Page = [0,'phs']
1412   
1413    G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.ConstraintMenu)
1414    SetStatusLine('')
1415   
1416    G2frame.Bind(wx.EVT_MENU, OnAddConstraint, id=G2G.wxID_CONSTRAINTADD)
1417    G2frame.Bind(wx.EVT_MENU, OnAddFunction, id=G2G.wxID_FUNCTADD)
1418    G2frame.Bind(wx.EVT_MENU, OnAddEquivalence, id=G2G.wxID_EQUIVADD)
1419    G2frame.Bind(wx.EVT_MENU, OnAddHold, id=G2G.wxID_HOLDADD)
1420    G2frame.Bind(wx.EVT_MENU, OnAddAtomEquiv, id=G2G.wxID_EQUIVALANCEATOMS)
1421#    G2frame.Bind(wx.EVT_MENU, OnAddRiding, id=G2G.wxID_ADDRIDING)
1422    def OnShowISODISTORT(event):
1423        ShowIsoDistortCalc(G2frame)
1424    G2frame.Bind(wx.EVT_MENU, OnShowISODISTORT, id=G2G.wxID_SHOWISO)
1425    # tab commands
1426    for id in (G2G.wxID_CONSPHASE,
1427               G2G.wxID_CONSHAP,
1428               G2G.wxID_CONSHIST,
1429               G2G.wxID_CONSGLOBAL,
1430               G2G.wxID_CONSSYM,
1431               ):
1432        G2frame.Bind(wx.EVT_MENU, RaisePage,id=id)
1433
1434    #G2frame.constr = G2G.GSNoteBook(parent=G2frame.dataWindow,size=G2frame.dataWindow.GetClientSize())
1435    G2frame.constr = G2G.GSNoteBook(parent=G2frame.dataWindow)
1436    G2frame.dataWindow.GetSizer().Add(G2frame.constr,1,wx.ALL|wx.EXPAND)
1437    # note that order of pages is hard-coded in RaisePage
1438    PhaseConstr = wx.ScrolledWindow(G2frame.constr)
1439    G2frame.constr.AddPage(PhaseConstr,'Phase')
1440    HAPConstr = wx.ScrolledWindow(G2frame.constr)
1441    G2frame.constr.AddPage(HAPConstr,'Histogram/Phase')
1442    HistConstr = wx.ScrolledWindow(G2frame.constr)
1443    G2frame.constr.AddPage(HistConstr,'Histogram')
1444    GlobalConstr = wx.ScrolledWindow(G2frame.constr)
1445    G2frame.constr.AddPage(GlobalConstr,'Global')
1446    SymConstr = wx.ScrolledWindow(G2frame.constr)
1447    G2frame.constr.AddPage(SymConstr,'Sym-Generated')
1448    wx.CallAfter(OnPageChanged,None)
1449    G2frame.constr.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, OnPageChanged)
1450    # validate all the constrants -- should not see any errors here normally
1451    errmsg,warnmsg = ConstraintsCheck(data)
1452    if errmsg:
1453        G2frame.ErrorDialog('Constraint Error',
1454                            'Error in constraints:\n'+errmsg+'\nCheck console output for more information',
1455                            parent=G2frame)
1456        print (errmsg)
1457        print (G2mv.VarRemapShow([],True))
1458    elif warnmsg:
1459        print ('Unexpected contraint warning:\n'+warnmsg)
1460       
1461################################################################################
1462# check scale & phase fractions, create constraint if needed
1463################################################################################
1464def CheckAllScalePhaseFractions(G2frame):
1465    '''Check if scale factor and all phase fractions are refined without a constraint
1466    for all used histograms, if so, offer the user a chance to create a constraint
1467    on the sum of phase fractions
1468    '''
1469    histograms, phases = G2frame.GetUsedHistogramsAndPhasesfromTree()
1470    for i,hist in enumerate(histograms):
1471        CheckScalePhaseFractions(G2frame,hist,histograms,phases)
1472       
1473def CheckScalePhaseFractions(G2frame,hist,histograms,phases):
1474    '''Check if scale factor and all phase fractions are refined without a constraint
1475    for histogram hist, if so, offer the user a chance to create a constraint
1476    on the sum of phase fractions
1477    '''
1478    if G2frame.testSeqRefineMode():
1479        histStr = '*'
1480    else:
1481        histStr = str(histograms[hist]['hId'])
1482    # Is this powder?
1483    if not hist.startswith('PWDR '): return
1484    # do this only if the scale factor is varied
1485    if not histograms[hist]['Sample Parameters']['Scale'][1]: return
1486    # are all phase fractions varied in all used histograms?
1487    phaseCount = 0
1488    for p in phases:
1489        if hist not in phases[p]['Histograms']: continue
1490        if phases[p]['Histograms'][hist]['Use'] and not phases[p]['Histograms'][hist]['Scale'][1]:
1491            return
1492        else:
1493            phaseCount += 1
1494   
1495    # all phase fractions and scale factor varied, now scan through constraints
1496    sub = G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Constraints') 
1497    Constraints = G2frame.GPXtree.GetItemPyData(sub)
1498    for c in Constraints.get('HAP',[]):
1499        if c[-1] != 'c': continue
1500        if not c[-3]: continue
1501        if len(c[:-3]) != phaseCount: continue
1502        # got a constraint equation with right number of terms, is it on phase fractions for
1503        # the correct histogram?
1504        if all([(i[1].name == 'Scale' and i[1].varname().split(':')[1] == histStr) for i in c[:-3]]):
1505            # got a constraint, this is OK
1506            return
1507    dlg = wx.MessageDialog(G2frame,'You are refining the scale factor and all phase fractions for histogram #'+
1508        histStr+'. This will produce an unstable refinement. '+
1509        'Do you want to constrain the sum of phase fractions?','Create constraint?',wx.OK|wx.CANCEL)
1510    if dlg.ShowModal() != wx.ID_OK:
1511        dlg.Destroy()
1512        return
1513    dlg.Destroy()
1514
1515    constr = []
1516    for p in phases:
1517        if hist not in phases[p]['Histograms']: continue
1518        if not phases[p]['Histograms'][hist]['Use']: continue
1519        constr += [[1.0,G2obj.G2VarObj(':'.join((str(phases[p]['pId']),histStr,'Scale')))]]
1520    constr += [1.0,None,'c']
1521    Constraints['HAP'] += [constr]
1522       
1523################################################################################
1524#### Make nuclear/magnetic phase transition constraints - called by OnTransform in G2phsGUI
1525################################################################################       
1526       
1527def TransConstraints(G2frame,oldPhase,newPhase,Trans,Vec,atCodes):
1528    '''Add constraints for new magnetic phase created via transformation of old
1529    nuclear one
1530    NB: A = [G11,G22,G33,2*G12,2*G13,2*G23]
1531    '''
1532   
1533    def SetUniqAj(pId,iA,SGData):
1534        SGLaue = SGData['SGLaue']
1535        SGUniq = SGData['SGUniq']
1536        if SGLaue in ['m3','m3m']:
1537            if iA in [0,1,2]:
1538                parm = '%d::%s'%(pId,'A0')
1539            else:
1540                parm = None
1541        elif SGLaue in ['4/m','4/mmm']:
1542            if iA in [0,1]:
1543                parm = '%d::%s'%(pId,'A0')
1544            elif iA == 2:
1545                parm = '%d::%s'%(pId,'A2')
1546            else:
1547                parm = None
1548        elif SGLaue in ['6/m','6/mmm','3m1', '31m', '3']:
1549            if iA in [0,1,3]:
1550                parm = '%d::%s'%(pId,'A0')
1551            elif iA == 2:
1552                parm = '%d::%s'%(pId,'A2')
1553            else:
1554                parm = None
1555        elif SGLaue in ['3R', '3mR']:
1556            if ia in [0,1,2]:
1557                parm = '%d::%s'%(pId,'A0')
1558            else:
1559                parm = '%d::%s'%(pId,'A3')
1560        elif SGLaue in ['mmm',]:
1561            if iA in [0,1,2]:
1562                parm = '%d::A%s'%(pId,iA)
1563            else:
1564                parm = None
1565        elif SGLaue == '2/m':
1566            if iA in [0,1,2]:
1567                parm = '%d::A%s'%(pId,iA)
1568            elif iA == 3 and SGUniq == 'c':
1569                parm = '%d::A%s'%(pId,iA)
1570            elif iA == 4 and SGUniq == 'b':
1571                parm = '%d::A%s'%(pId,iA)
1572            elif iA == 5 and SGUniq == 'a':
1573                parm = '%d::A%s'%(pId,iA)
1574            else:
1575                parm = None           
1576        else:
1577            parm = '%d::A%s'%(pId,iA)
1578        return parm
1579   
1580    Histograms,Phases = G2frame.GetUsedHistogramsAndPhasesfromTree()
1581    UseList = newPhase['Histograms']
1582    detTrans = np.abs(nl.det(Trans))
1583    invTrans = nl.inv(Trans)
1584    opId = oldPhase['pId']
1585    npId = newPhase['pId']
1586    cx,ct,cs,cia = newPhase['General']['AtomPtrs']
1587    nAtoms = newPhase['Atoms']
1588    oSGData = oldPhase['General']['SGData']
1589    nSGData = newPhase['General']['SGData']
1590    #oAcof = G2lat.cell2A(oldPhase['General']['Cell'][1:7])
1591    #nAcof = G2lat.cell2A(newPhase['General']['Cell'][1:7])
1592    item = G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Constraints')
1593    if not item:
1594        print('Error: no constraints in Data Tree')
1595        return
1596    constraints = G2frame.GPXtree.GetItemPyData(item)
1597    xnames = ['dAx','dAy','dAz']
1598    # constraints on matching atom params between phases
1599    for ia,code in enumerate(atCodes):
1600        atom = nAtoms[ia]
1601        if not ia and atom[cia] == 'A':
1602            wx.MessageDialog(G2frame,
1603                'Anisotropic thermal motion constraints are not developed at the present time',
1604                'Anisotropic thermal constraint?',style=wx.ICON_INFORMATION).ShowModal()
1605        siteSym = G2spc.SytSym(atom[cx:cx+3],nSGData)[0]
1606        CSX = G2spc.GetCSxinel(siteSym)
1607#        CSU = G2spc.GetCSuinel(siteSym)
1608        item = code.split('+')[0]
1609        iat,opr = item.split(':')
1610        Nop = abs(int(opr))%100-1
1611        if '-' in opr:
1612            Nop *= -1
1613        Opr = oldPhase['General']['SGData']['SGOps'][abs(Nop)][0]
1614        if Nop < 0:         #inversion
1615            Opr *= -1
1616        XOpr = np.inner(Opr,Trans)
1617        for i,ix in enumerate(list(CSX[0])):
1618            if not ix:
1619                continue
1620            name = xnames[i]
1621            IndpCon = [1.0,G2obj.G2VarObj('%d::%s:%d'%(npId,name,ia))]
1622            DepCons = []
1623            for iop,opval in enumerate(XOpr[i]):
1624                if opval:
1625                    DepCons.append([opval,G2obj.G2VarObj('%d::%s:%s'%(opId,xnames[iop],iat))])
1626            if len(DepCons) == 1:
1627                constraints['Phase'].append([DepCons[0],IndpCon,None,None,'e'])
1628            elif len(DepCons) > 1:
1629                for Dep in DepCons:
1630                    Dep[0] *= -1
1631                constraints['Phase'].append([IndpCon]+DepCons+[0.0,None,'c'])
1632        for name in ['Afrac','AUiso']:
1633            IndpCon = [1.0,G2obj.G2VarObj('%d::%s:%d'%(npId,name,ia))]
1634            DepCons = [1.0,G2obj.G2VarObj('%d::%s:%s'%(opId,name,iat))]
1635            constraints['Phase'].append([DepCons,IndpCon,None,None,'e'])
1636           
1637        # unfinished Anisotropic constraint generation
1638#        Uids = [[0,0,'AU11'],[1,1,'AU22'],[2,2,'AU33'],[0,1,'AU12'],[0,2,'AU13'],[1,2,'AU23']]
1639#        DepConsDict = dict(zip(Us,[[],[],[],[],[],[]]))
1640#        for iu,Uid in enumerate(Uids):
1641#            UMT = np.zeros((3,3))
1642#            UMT[Uid[0],Uid[1]] = 1
1643#            nUMT = G2lat.prodMGMT(UMT,invTrans)
1644#            nUT = G2lat.UijtoU6(nUMT)
1645#            for iu,nU in enumerate(nUT):
1646#                if abs(nU) > 1.e-8:
1647#                    parm = '%d::%s;%s'%(opId,Us[iu],iat)
1648#                    DepConsDict[Uid[2]].append([abs(nU%1.),G2obj.G2VarObj(parm)])
1649#        nUcof = atom[iu:iu+6]
1650#        conStrings = []
1651#        for iU,Usi in enumerate(Us):
1652#            parm = '%d::%s;%d'%(npId,Usi,ia)
1653#            parmDict[parm] = nUcof[iU]
1654#            varyList.append(parm)
1655#            IndpCon = [1.0,G2obj.G2VarObj(parm)]
1656#            conStr = str([IndpCon,DepConsDict[Usi]])
1657#            if conStr in conStrings:
1658#                continue
1659#            conStrings.append(conStr)
1660#            if len(DepConsDict[Usi]) == 1:
1661#                if DepConsDict[Usi][0]:
1662#                    constraints['Phase'].append([IndpCon,DepConsDict[Usi][0],None,None,'e'])
1663#            elif len(DepConsDict[Usi]) > 1:       
1664#                for Dep in DepConsDict[Usi]:
1665#                    Dep[0] *= -1
1666#                constraints['Phase'].append([IndpCon]+DepConsDict[Usi]+[0.0,None,'c'])
1667           
1668        #how do I do Uij's for most Trans?
1669
1670    # constraints on lattice parameters between phases
1671#    T = nl.inv(Trans).T
1672    T = Trans.T
1673    conMat = [
1674        [T[0,0]**2,T[0,1]**2,T[0,2]**2,T[0,0]*T[0,1],T[0,0]*T[0,2],T[0,1]*T[0,2]],
1675        [T[1,0]**2,T[1,1]**2,T[1,2]**2,T[1,0]*T[1,1],T[1,0]*T[1,2],T[1,1]*T[1,2]],
1676        [T[2,0]**2,T[2,1]**2,T[2,2]**2,T[2,0]*T[2,1],T[2,0]*T[2,2],T[2,1]*T[2,2]],
1677        [2.*T[0,0]*T[1,0],2.*T[0,1]*T[1,1],2.*T[0,2]*T[1,2],T[0,0]*T[1,1]+T[0,1]*T[1,0],T[0,0]*T[1,2]+T[0,2]*T[1,0],T[0,1]*T[1,2]+T[0,2]*T[1,1]],
1678        [2.*T[0,0]*T[2,0],2.*T[0,1]*T[2,1],2.*T[0,2]*T[2,2],T[0,0]*T[2,1]+T[0,1]*T[2,0],T[0,0]*T[2,2]+T[0,2]*T[2,0],T[0,1]*T[2,2]+T[0,2]*T[2,1]],
1679        [2.*T[1,0]*T[2,0],2.*T[1,1]*T[2,1],2.*T[1,2]*T[2,2],T[1,0]*T[2,1]+T[1,1]*T[2,0],T[1,0]*T[2,2]+T[1,2]*T[2,0],T[1,1]*T[2,2]+T[1,2]*T[2,1]]
1680        ]
1681    # Gnew = conMat * A:
1682#         T00**2*a0  T01**2*a1 T02**2*a2 T00*T01*a3    T00*T02*a4    T01*T02*a5
1683#         T10**2*a0  T11**2*a1 T12**2*a2 T10*T11*a3    T10*T12*a4    T11*T12*a5
1684#         T20**2*a0  T21**2*a1 T22**2*a2 T20*T21*a3    T20*T22*a4    T21*T22*a5
1685#         2*T00*T10*a0      2*T01*T11*a1     2*T02*T12*a2     (T00*T11 + T01*T10)*a3      (T00*T12 + T02*T10)*a4      (T01*T12 + T02*T11)*a5
1686#         2*T00*T20*a0      2*T01*T21*a1     2*T02*T22*a2     (T00*T21 + T01*T20)*a3      (T00*T22 + T02*T20)*a4      (T01*T22 + T02*T21)*a5
1687#         2*T10*T20*a0      2*T11*T21*a1     2*T12*T22*a2     (T10*T21 + T11*T20)*a3      (T10*T22 + T12*T20)*a4      (T11*T22 + T12*T21)*a5
1688    # Generated as symbolic code using:
1689    # import sym
1690    # A0, A1, A2, A3, A4, A5 = sym.symbols('A0, A1, A2, A3, A4, A5')
1691    # G = sym.Matrix([ [A0,    A3/2.,  A4/2.], [A3/2.,  A1,    A5/2.], [A4/2., A5/2.,    A2]])
1692    # transformation matrix
1693    # T00, T10, T20, T01, T11, T21, T02, T12, T22 = sym.symbols('T00, T10, T20, T01, T11, T21, T02, T12, T22')
1694    # Tr = sym.Matrix([ [T00, T10, T20], [T01, T11, T21], [T02, T12, T22],])
1695    # Gnew = (Tr.T*G)*Tr
1696   
1697    #print('old A',G2lat.cell2A(oldPhase['General']['Cell'][1:7]))
1698    #print('new A',G2lat.cell2A(newPhase['General']['Cell'][1:7]))
1699   
1700#this is still incorrect for hex/trig/ortho/tetragonal --> monoclinic
1701   
1702#    for iAnew,Asi in enumerate(['A0','A1','A2','A3','A4','A5']): # loop through A[i] for new cell
1703#        Nparm = str(npId) + '::' + Asi
1704#        if Nparm != SetUniqAj(npId,iAnew,nSGData):
1705#            continue # skip: Ai constrained from Aj or must be zero
1706#        multDict = {}
1707#        for iAorg in range(6):
1708#            cA = conMat[iAnew][iAorg] # coeff for A[i] in constraint matrix
1709#            if abs(cA) < 1.e-8: continue
1710#            parm = SetUniqAj(opId,iAorg,oSGData) # translate to unique A[i] in original cell
1711#            if not parm: continue # must be zero
1712#            # sum coeff
1713#            if parm in multDict:
1714#                multDict[parm] += cA
1715#            else:
1716#                multDict[parm] = cA
1717#        # any non-zero multipliers?
1718#        maxMult = 0
1719#        for i in multDict:
1720#            maxMult = max(maxMult,abs(multDict[i]))
1721#        if maxMult <= 0:  # Nparm computes as zero; Fix this parameter
1722#            constraints['Phase'] += [[
1723#                [0.0,G2obj.G2VarObj(Nparm)],
1724#                None,None,'h']]
1725#        elif len(multDict) == 1:        # create equivalence
1726#            key = list(multDict.keys())[0]
1727#            constraints['Phase'] += [[
1728#                [1.0,G2obj.G2VarObj(key)],
1729#                [multDict[key],G2obj.G2VarObj(Nparm)],
1730#                None,None,'e']]
1731#        else:                           # create constraint
1732#            constr = [[-1.0,G2obj.G2VarObj(Nparm)]]
1733#            for key in multDict:
1734#                constr += [[multDict[key],G2obj.G2VarObj(key)]]
1735#            constr += [0.0,None,'c']
1736#            constraints['Phase'] += [constr]
1737   
1738    # constraints on HAP Scale, etc.
1739    for hId,hist in enumerate(UseList):    #HAP - seems OK
1740        ohapkey = '%d:%d:'%(opId,hId)
1741        nhapkey = '%d:%d:'%(npId,hId)
1742        IndpCon = [1.0,G2obj.G2VarObj(ohapkey+'Scale')]
1743        DepCons = [detTrans,G2obj.G2VarObj(nhapkey+'Scale')]
1744        constraints['HAP'].append([DepCons,IndpCon,None,None,'e'])
1745        for name in ['Size;i','Mustrain;i']:
1746            IndpCon = [1.0,G2obj.G2VarObj(ohapkey+name)]
1747            DepCons = [1.0,G2obj.G2VarObj(nhapkey+name)]
1748            constraints['HAP'].append([IndpCon,DepCons,None,None,'e'])
1749       
1750################################################################################
1751#### Rigid bodies
1752################################################################################
1753resRBsel = None
1754def UpdateRigidBodies(G2frame,data):
1755    '''Called when Rigid bodies tree item is selected.
1756    Displays the rigid bodies in the data window
1757    '''
1758    if not data.get('RBIds') or not data:
1759        data.update({'Vector':{'AtInfo':{}},'Residue':{'AtInfo':{}},
1760            'RBIds':{'Vector':[],'Residue':[]}})       #empty/bad dict - fill it
1761           
1762    global resList,resRBsel
1763    Indx = {}
1764    resList = []
1765    plotDefaults = {'oldxy':[0.,0.],'Quaternion':[0.,0.,0.,1.],'cameraPos':30.,'viewDir':[0,0,1],}
1766
1767    G2frame.rbBook = G2G.GSNoteBook(parent=G2frame.dataWindow)
1768    G2frame.dataWindow.GetSizer().Add(G2frame.rbBook,1,wx.ALL|wx.EXPAND)
1769    VectorRB = wx.ScrolledWindow(G2frame.rbBook)
1770    VectorRBDisplay = wx.Panel(VectorRB)
1771    G2frame.rbBook.AddPage(VectorRB,'Vector rigid bodies')
1772    ResidueRB = wx.ScrolledWindow(G2frame.rbBook)
1773    ResidueRBDisplay = wx.Panel(ResidueRB)
1774    G2frame.rbBook.AddPage(ResidueRB,'Residue rigid bodies')
1775   
1776    def OnPageChanged(event):
1777        global resList
1778        resList = []
1779        if event:       #page change event!
1780            page = event.GetSelection()
1781        else:
1782            page = G2frame.rbBook.GetSelection()
1783        G2frame.rbBook.ChangeSelection(page)
1784        text = G2frame.rbBook.GetPageText(page)
1785        if text == 'Vector rigid bodies':
1786            G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.VectorBodyMenu)
1787            G2frame.Bind(wx.EVT_MENU, AddVectorRB, id=G2G.wxID_VECTORBODYADD)
1788            G2frame.Bind(wx.EVT_MENU, ExtractPhaseRB, id=G2G.wxID_VECTORBODYIMP)
1789            G2frame.Bind(wx.EVT_MENU, AddVectTrans, id=G2G.wxID_VECTORBODYEXTD)
1790            G2frame.Bind(wx.EVT_MENU, SaveVectorRB, id=G2G.wxID_VECTORBODYSAV)
1791            G2frame.Bind(wx.EVT_MENU, ReadVectorRB, id=G2G.wxID_VECTORBODYRD)
1792            G2frame.Page = [page,'vrb']
1793            UpdateVectorRB()
1794        elif text == 'Residue rigid bodies':
1795            G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.RigidBodyMenu)
1796            G2frame.Bind(wx.EVT_MENU, AddResidueRB, id=G2G.wxID_RIGIDBODYADD)
1797            G2frame.Bind(wx.EVT_MENU, ExtractPhaseRB, id=G2G.wxID_RIGIDBODYIMP)
1798            G2frame.Bind(wx.EVT_MENU, OnImportRigidBody, id=G2G.wxID_RIGIDBODYIMPORT)
1799            G2frame.Bind(wx.EVT_MENU, OnSaveRigidBody, id=G2G.wxID_RIGIDBODYSAVE)
1800            G2frame.Bind(wx.EVT_MENU, OnDefineTorsSeq, id=G2G.wxID_RESIDUETORSSEQ) #enable only if residue RBs exist?
1801            G2frame.Bind(wx.EVT_MENU, DumpVectorRB, id=G2G.wxID_RESBODYSAV)
1802            G2frame.Bind(wx.EVT_MENU, LoadVectorRB, id=G2G.wxID_RESBODYRD)
1803            G2frame.Page = [page,'rrb']
1804            UpdateResidueRB()
1805        else:
1806            G2gd.SetDataMenuBar(G2frame)
1807            #G2frame.Page = [page,'rrb']
1808           
1809    def getMacroFile(macName):
1810        defDir = os.path.join(os.path.split(__file__)[0],'GSASIImacros')
1811        dlg = wx.FileDialog(G2frame,message='Choose '+macName+' rigid body macro file',
1812            defaultDir=defDir,defaultFile="",wildcard="GSAS-II macro file (*.mac)|*.mac",
1813            style=wx.FD_OPEN | wx.FD_CHANGE_DIR)
1814        try:
1815            if dlg.ShowModal() == wx.ID_OK:
1816                macfile = dlg.GetPath()
1817                macro = open(macfile,'Ur')
1818                head = macro.readline()
1819                if macName not in head:
1820                    print (head)
1821                    print ('**** ERROR - wrong restraint macro file selected, try again ****')
1822                    macro = []
1823            else: # cancel was pressed
1824                macro = []
1825        finally:
1826            dlg.Destroy()
1827        return macro        #advanced past 1st line
1828       
1829    def getTextFile():
1830        dlg = wx.FileDialog(G2frame,'Choose rigid body text file', G2frame.LastGPXdir, '',
1831            "GSAS-II text file (*.txt)|*.txt|XYZ file (*.xyz)|*.xyz|"
1832            "Sybyl mol2 file (*.mol2)|*.mol2|PDB file (*.pdb;*.ent)|*.pdb;*.ent",
1833            wx.FD_OPEN | wx.FD_CHANGE_DIR)
1834        try:
1835            if dlg.ShowModal() == wx.ID_OK:
1836                txtfile = dlg.GetPath()
1837                ext = os.path.splitext(txtfile)[1]
1838                text = open(txtfile,'Ur')
1839            else: # cancel was pressed
1840                ext = ''
1841                text = []
1842        finally:
1843            dlg.Destroy()
1844        if 'ent' in ext:
1845            ext = '.pdb'
1846        return text,ext.lower()
1847       
1848    def OnImportRigidBody(event):
1849        page = G2frame.rbBook.GetSelection()
1850        if 'Vector' in G2frame.rbBook.GetPageText(page):
1851            pass
1852        elif 'Residue' in G2frame.rbBook.GetPageText(page):
1853            ImportResidueRB()
1854           
1855    def OnSaveRigidBody(event):
1856        page = G2frame.rbBook.GetSelection()
1857        if 'Vector' in G2frame.rbBook.GetPageText(page):
1858            pass
1859        elif 'Residue' in G2frame.rbBook.GetPageText(page):
1860            SaveResidueRB()
1861           
1862    def DumpVectorRB(event):
1863        global resRBsel
1864        if resRBsel not in data['Residue']:
1865            return
1866        rbData = data['Residue'][resRBsel]
1867        pth = G2G.GetExportPath(G2frame)
1868        dlg = wx.FileDialog(G2frame, 'Choose file to save residue rigid body',
1869            pth, '', 'RRB files (*.resbody)|*.resbody',
1870            wx.FD_SAVE|wx.FD_OVERWRITE_PROMPT)
1871        try:
1872            if dlg.ShowModal() == wx.ID_OK:
1873                filename = dlg.GetPath()
1874                filename = os.path.splitext(filename)[0]+'.resbody'  # set extension
1875                fp = open(filename,'w')
1876                fp.write('Name: '+rbData['RBname']+'\n')
1877                fp.write('atNames: ')
1878                for i in rbData['atNames']:
1879                    fp.write(str(i)+" ") 
1880                fp.write('\n')
1881                for item in rbData['rbSeq']:
1882                    fp.write('rbSeq: ') 
1883                    fp.write('{:d} {:d} {:.1f}: '.format(*item[:3]))
1884                    for num in item[3]:
1885                        fp.write('{:d} '.format(num))
1886                    fp.write('\n')
1887                for i,sym in enumerate(rbData['rbTypes']):
1888                    fp.write("{:3s}".format(sym))
1889                    fp.write('{:8.5f}{:9.5f}{:9.5f}   '
1890                            .format(*rbData['rbXYZ'][i]))
1891                    fp.write('\n')
1892                fp.close()
1893                print ('Vector rigid body saved to: '+filename)
1894               
1895        finally:
1896            dlg.Destroy()
1897       
1898    def LoadVectorRB(event):
1899        AtInfo = data['Residue']['AtInfo']
1900        pth = G2G.GetExportPath(G2frame)
1901        dlg = wx.FileDialog(G2frame, 'Choose file to read vector rigid body',
1902            pth, '', 'RRB files (*.resbody)|*.resbody',
1903            wx.FD_OPEN)
1904        try:
1905            if dlg.ShowModal() == wx.ID_OK:
1906                filename = dlg.GetPath()
1907                filename = os.path.splitext(filename)[0]+'.resbody'  # set extension
1908                fp = open(filename,'r')
1909                l = fp.readline().strip()
1910                if 'Name' not in l:
1911                    fp.close()
1912                    G2frame.ErrorDialog('Read Error',
1913                        'File '+filename+' does not start with Name\nFirst line ='
1914                        +l+'\ninvalid file',parent=G2frame)
1915                    return
1916                name = l.split(':')[1].strip()
1917                line = fp.readline().strip().split(':')[1].split()
1918                atNames = [i for i in line]
1919                types = []
1920                coords = []
1921                l = fp.readline().strip()
1922                rbSeq = []
1923                while 'rbSeq' in l:
1924                    tag,vals,lst = l.split(':')
1925                    seq = []
1926                    for t,v in zip((int,int,float),vals.split()): 
1927                        seq.append(t(v))
1928                    seq.append([])
1929                    for num in lst.split():
1930                        seq[-1].append(int(num))
1931                    rbSeq.append(seq)
1932                    l = fp.readline().strip()                   
1933                while l:
1934                    nums = l.strip().split()
1935                    types.append(nums.pop(0))
1936                    t = types[-1]
1937                    if t not in AtInfo:
1938                        Info = G2elem.GetAtomInfo(t)
1939                        AtInfo[t] = [Info['Drad'],Info['Color']]
1940                    coords.append([float(nums.pop(0)) for j in range(3)])
1941                    l = fp.readline().strip()
1942                fp.close()
1943            else:
1944                return       
1945        finally:
1946            dlg.Destroy()
1947        coords = np.array(coords)
1948        rbid = ran.randint(0,sys.maxsize)
1949        data['Residue'][rbid] = {'RBname':name,
1950                'rbXYZ': coords,
1951                'rbRef':[0,1,2,False],
1952                'rbTypes':types, 'atNames':atNames,
1953                'useCount':0,
1954                'rbSeq':rbSeq, 'SelSeq':[0,0],}
1955        data['RBIds']['Residue'].append(rbid)
1956        UpdateResidueRB()
1957   
1958    def AddVectorRB(event):
1959        'Create a new vector rigid body'
1960        AtInfo = data['Vector']['AtInfo']
1961        dlg = G2G.MultiIntegerDialog(G2frame,'New Rigid Body',['No. atoms','No. translations'],[1,1])
1962        if dlg.ShowModal() == wx.ID_OK:
1963            nAtoms,nTrans = dlg.GetValues()
1964            rbid = ran.randint(0,sys.maxsize)
1965            vecMag = [1.0 for i in range(nTrans)]
1966            vecRef = [False for i in range(nTrans)]
1967            vecVal = [np.zeros((nAtoms,3)) for j in range(nTrans)]
1968            rbTypes = ['C' for i in range(nAtoms)]
1969            Info = G2elem.GetAtomInfo('C')
1970            AtInfo['C'] = [Info['Drad'],Info['Color']]
1971            data['Vector'][rbid] = {'RBname':'UNKRB','VectMag':vecMag,'rbXYZ':np.zeros((nAtoms,3)),
1972                'rbRef':[0,1,2,False],'VectRef':vecRef,'rbTypes':rbTypes,'rbVect':vecVal,'useCount':0}
1973            data['RBIds']['Vector'].append(rbid)
1974        dlg.Destroy()
1975        UpdateVectorRB()
1976
1977    def ExtractPhaseRB(event):
1978        'Extract a rigid body from a file with a phase'
1979        def SetupDrawing(atmData):
1980            '''Add the dicts needed for G2plt.PlotStructure to work to the
1981            reader .Phase object
1982            '''
1983            generalData = atmData['General']
1984            generalData['BondRadii'] = []
1985
1986            G2phG.SetDrawingDefaults(atmData['Drawing'])
1987            atmData['Drawing'].update(
1988                {'oldxy':[0.,0.],'Quaternion':[0.,0.,0.,1.],'cameraPos':150.,
1989                     'viewDir':[0,0,1],'atomPtrs': [2, 1, 6, 17],
1990                     })
1991            atmData['Drawing']['showRigidBodies'] = False
1992            generalData['Map'] = {'MapType':False, 'rho':[]}
1993            generalData['AtomTypes'] = []
1994            generalData['BondRadii'] = []
1995            generalData['AngleRadii'] = []
1996            generalData['vdWRadii'] = []
1997            generalData['Color'] = []
1998            generalData['Isotopes'] = {}
1999            generalData['Isotope'] = {}
2000            cx,ct,cs,cia = generalData['AtomPtrs']
2001            generalData['Mydir'] = G2frame.dirname
2002            for iat,atom in enumerate(atmData['Atoms']):
2003                atom[ct] = atom[ct].lower().capitalize()      #force elem symbol to standard form
2004                if atom[ct] not in generalData['AtomTypes'] and atom[ct] != 'UNK':
2005                    Info = G2elem.GetAtomInfo(atom[ct])
2006                    if not Info:
2007                        atom[ct] = 'UNK'
2008                        continue
2009                    atom[ct] = Info['Symbol'] # N.B. symbol might be changed by GetAtomInfo
2010                    generalData['AtomTypes'].append(atom[ct])
2011                    generalData['Z'] = Info['Z']
2012                    generalData['Isotopes'][atom[ct]] = Info['Isotopes']
2013                    generalData['BondRadii'].append(Info['Drad'])
2014                    generalData['AngleRadii'].append(Info['Arad'])
2015                    generalData['vdWRadii'].append(Info['Vdrad'])
2016                    if atom[ct] in generalData['Isotope']:
2017                        if generalData['Isotope'][atom[ct]] not in generalData['Isotopes'][atom[ct]]:
2018                            isotope = list(generalData['Isotopes'][atom[ct]].keys())[-1]
2019                            generalData['Isotope'][atom[ct]] = isotope
2020                    else:
2021                        generalData['Isotope'][atom[ct]] = 'Nat. Abund.'
2022                        if 'Nat. Abund.' not in generalData['Isotopes'][atom[ct]]:
2023                            isotope = list(generalData['Isotopes'][atom[ct]].keys())[-1]
2024                            generalData['Isotope'][atom[ct]] = isotope
2025                    generalData['Color'].append(Info['Color'])
2026                    if generalData['Type'] == 'magnetic':
2027                        if len(landeg) < len(generalData['AtomTypes']):
2028                            landeg.append(2.0)
2029            atmData['Drawing']['Atoms'] = []
2030            for atom in atmData['Atoms']:
2031                atmData['Drawing']['Atoms'].append(MakeDrawAtom(atmData,atom))
2032
2033        def onCancel(event,page=0):
2034            'complete or bail out from RB define, cleaning up'
2035            G2frame.rbBook.DeletePage(G2frame.rbBook.FindPage(pagename))
2036            G2frame.rbBook.SetSelection(page)
2037
2038        def Page1():
2039            '''Show the GUI for first stage of the rigid body with all atoms in
2040            phase in crystal coordinates. Select the atoms to go onto the
2041            next stage
2042            '''
2043            def ShowSelection(selections):
2044                'respond to change in atom selections'
2045                ct,cs = [1,8]
2046                generalData = rd.Phase['General']
2047                for i,atom in enumerate(rd.Phase['Drawing']['Atoms']):
2048                    if i in selections:
2049                        factor = 1
2050                    else:
2051                        factor = 2.5
2052                    atNum = generalData['AtomTypes'].index(atom[ct]) 
2053                    atom[cs] = list(np.array(generalData['Color'][atNum])//factor) 
2054                draw(*drawArgs)
2055            def onPage1OK(event):
2056                '1st section has been completed, move onto next'
2057                G2frame.G2plotNB.Delete(rd.Phase['General']['Name'])
2058                GetCoords(atmsel)
2059                Page2()
2060
2061            SetupDrawing(rd.Phase) # add information to reader object to allow plotting
2062            atomlist = [atom[0] for atom in rd.Phase['Atoms']]
2063            atmsel = list(range(len(rd.Phase['Atoms'])))
2064            # broken -- # why no bonds?
2065            #for atm in rd.Phase['Drawing']['Atoms']:
2066            #    atm[6] = 'balls & sticks'
2067
2068            draw,drawArgs = G2plt.PlotStructure(G2frame,rd.Phase,True)
2069            ShowSelection(atmsel)
2070
2071            if G2frame.rbBook.FindPage(pagename) is not None:
2072                G2frame.rbBook.DeletePage(G2frame.rbBook.FindPage(pagename))
2073
2074            RBImp = wx.ScrolledWindow(G2frame.rbBook)
2075            RBImpPnl = wx.Panel(RBImp)
2076            G2frame.rbBook.AddPage(RBImp,pagename)
2077            G2frame.rbBook.SetSelection(G2frame.rbBook.FindPage(pagename))
2078
2079            mainSizer = G2G.G2MultiChoiceWindow(RBImpPnl,
2080                            'Select atoms to import',
2081                            atomlist,atmsel,OnChange=ShowSelection)
2082
2083            # OK/Cancel buttons       
2084            btnsizer = wx.StdDialogButtonSizer()
2085            OKbtn = wx.Button(RBImpPnl, wx.ID_OK, 'Continue')
2086            OKbtn.SetDefault()
2087            btnsizer.AddButton(OKbtn)
2088            OKbtn.Bind(wx.EVT_BUTTON,onPage1OK)
2089            btn = wx.Button(RBImpPnl, wx.ID_CANCEL)
2090            btn.Bind(wx.EVT_BUTTON,onCancel)
2091            btnsizer.AddButton(btn)
2092            btnsizer.Realize()
2093            mainSizer.Add(btnsizer,0,wx.ALIGN_CENTER,50)
2094
2095            RBImpPnl.SetSizer(mainSizer,True)
2096
2097            mainSizer.Layout()   
2098            Size = mainSizer.GetMinSize()
2099            Size[0] += 40
2100            Size[1] = max(Size[1],G2frame.GetSize()[1]-200) + 20
2101            RBImpPnl.SetSize(Size)
2102            RBImp.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
2103            RBImp.Scroll(0,0)
2104
2105        def Page2():
2106            '''Show the GUI for the second stage, where selected atoms are
2107            now in Cartesian space, manipulate the axes and export selected
2108            atoms to a vector or residue rigid body.
2109            '''
2110            def UpdateDraw(event=None):
2111                'Called when info changes in grid, replots'
2112                UpdateVectorBody(rbData)
2113                DrawCallback()
2114               
2115            def onSetAll(event):
2116                'Set all atoms as selected'
2117                grid.completeEdits()
2118                for i in range(len(rd.Phase['RBselection'])):
2119                    rd.Phase['RBselection'][i] = 1 # table needs 0/1 for T/F
2120                grid.ForceRefresh()
2121                UpdateDraw()
2122               
2123            def onToggle(event):
2124                'Toggles selection state for all atoms'
2125                grid.completeEdits()
2126                for i in range(len(rd.Phase['RBselection'])):
2127                    rd.Phase['RBselection'][i] = int(not rd.Phase['RBselection'][i])
2128                grid.ForceRefresh()
2129                UpdateDraw()
2130               
2131            def onSetOrigin(event):
2132                'Resets origin to midpoint between all selected atoms'
2133                grid.completeEdits()
2134                center = np.array([0.,0.,0.])
2135                count = 0
2136                for i in range(len(rd.Phase['RBselection'])):
2137                    if rd.Phase['RBselection'][i]:
2138                        count += 1
2139                        center += rd.Phase['RBcoords'][i]
2140                if count:
2141                    rd.Phase['RBcoords'] -= center/count
2142                grid.ForceRefresh()
2143                UpdateDraw()
2144               
2145            def onSetX(event):
2146                grid.completeEdits()
2147                center = np.array([0.,0.,0.])
2148                count = 0
2149                for i in range(len(rd.Phase['RBselection'])):
2150                    if rd.Phase['RBselection'][i]:
2151                        count += 1
2152                        center += rd.Phase['RBcoords'][i]
2153                if not count:
2154                    G2G.G2MessageBox(G2frame,'No atoms selected',
2155                                    'Selection required')
2156                    return
2157                XYZP = center/count
2158                if np.sqrt(sum(XYZP**2)) < 0.1:
2159                    G2G.G2MessageBox(G2frame,
2160                            'The selected atom(s) are too close to the origin',
2161                            'near origin')
2162                    return
2163                if bntOpts['direction'] == 'y':
2164                    YP = XYZP / np.sqrt(np.sum(XYZP**2))
2165                    ZP = np.cross((1,0,0),YP)
2166                    if sum(ZP*ZP) < .1: # pathological condition: Y' along X
2167                        ZP = np.cross((0,0,1),YP)
2168                    XP = np.cross(YP,ZP)
2169                elif bntOpts['direction'] == 'z':
2170                    ZP = XYZP / np.sqrt(np.sum(XYZP**2))
2171                    XP = np.cross((0,1,0),ZP)
2172                    if sum(XP*XP) < .1: # pathological condition: X' along Y
2173                        XP = np.cross((0,0,1),ZP)
2174                    YP = np.cross(ZP,XP)
2175                else:
2176                    XP = XYZP / np.sqrt(np.sum(XYZP**2))
2177                    YP = np.cross((0,0,1),XP)
2178                    if sum(YP*YP) < .1: # pathological condition: X' along Z
2179                        YP = np.cross((0,1,0),XP)
2180                    ZP = np.cross(XP,YP)
2181                trans = np.array((XP,YP,ZP))
2182                # update atoms in place
2183                rd.Phase['RBcoords'][:] = np.inner(trans,rd.Phase['RBcoords']).T
2184                grid.ForceRefresh()
2185                UpdateDraw()
2186
2187            def onSetPlane(event): 
2188                '''Compute least-squares plane for selected atoms;
2189                move atoms so that LS plane aligned with x-y plane,
2190                with minimum required change to x
2191                '''
2192                grid.completeEdits()
2193                selList = [i==1 for i in rd.Phase['RBselection']]
2194                XYZ = rd.Phase['RBcoords'][selList]
2195                if len(XYZ) < 3: 
2196                    G2G.G2MessageBox(G2frame,'A plane requires three or more atoms',
2197                                     'Need more atoms')
2198                    return
2199                # fit 3 ways (in case of singularity) and take result with lowest residual
2200                X,Y,Z = [XYZ[:,i] for i in (0,1,2)]
2201                XZ = copy.copy(XYZ)
2202                XZ[:,1] = 1
2203                (a,d,b), resd, rank, sing = nl.lstsq(XZ, -Y)
2204                resid_min = resd
2205                normal = a,1,b
2206                YZ = copy.copy(XYZ)
2207                YZ[:,0] = 1
2208                (d,a,b), resd, rank, sing = nl.lstsq(YZ, -X)
2209                if resid_min > resd:
2210                    resid_min = resd
2211                    normal = 1,a,b
2212                XY = copy.copy(XYZ)
2213                XY[:,2] = 1
2214                (a,b,d), resd, rank, sing = nl.lstsq(XY, -Z)
2215                if resid_min > resd:
2216                    resid_min = resd
2217                    normal = a,b,1
2218                # solve for  ax + bx + z + c = 0 or equivalently ax + bx + c = -z
2219                # try:
2220                # except:
2221                #     G2G.G2MessageBox(G2frame,
2222                #             'Error computing plane; are atoms in a line?',
2223                #             'Computation error')
2224                #     return
2225                if bntOpts['plane'] == 'xy':
2226                    # new coordinate system is
2227                    #   ZP, z' normal to plane
2228                    #   YP, y' = z' cross x (= [0,1,-b])
2229                    #   XP, x' = (z' cross x) cross z'
2230                    # this puts XP as close as possible to X with XP & YP in plane
2231                    ZP = np.array(normal)
2232                    ZP /= np.sqrt(np.sum(ZP**2))
2233                    YP = np.cross(ZP,[1,0,0])
2234                    if sum(YP*YP) < .1: # pathological condition: z' along x
2235                        YP = np.cross(ZP,[0,1,0])
2236                    YP /= np.sqrt(np.sum(YP**2))
2237                    XP = np.cross(YP,ZP)
2238                elif bntOpts['plane'] == 'yz':
2239                    # new coordinate system is
2240                    #   XP, x' normal to plane
2241                    #   ZP, z' = x' cross y
2242                    #   YP, y' = (x' cross y) cross x'
2243                    # this puts y' as close as possible to y with z' & y' in plane
2244                    XP = np.array(normal)
2245                    XP /= np.sqrt(np.sum(XP**2))
2246                    ZP = np.cross(XP,[0,1,0])
2247                    if sum(ZP*ZP) < .1: # pathological condition: x' along y
2248                        ZP = np.cross(XP,(0,0,1))
2249                    ZP /= np.sqrt(np.sum(ZP**2))
2250                    YP = np.cross(ZP,XP)
2251                elif bntOpts['plane'] == 'xz':
2252                    # new coordinate system is
2253                    #   YP, y' normal to plane
2254                    #   ZP, z' = x cross y'
2255                    #   XP, y' = (x cross y') cross z'
2256                    # this puts XP as close as possible to X with XP & YP in plane
2257                    YP = np.array(normal)
2258                    YP /= np.sqrt(np.sum(YP**2))
2259                    ZP = np.cross([1,0,0],YP)
2260                    if sum(ZP*ZP) < .1: # pathological condition: y' along x
2261                        ZP = np.cross([0,1,0],YP)
2262                    ZP /= np.sqrt(np.sum(ZP**2))
2263                    XP = np.cross(YP,ZP)
2264                else:
2265                    print('unexpected plane',bntOpts['plane'])
2266                    return
2267                trans = np.array((XP,YP,ZP))
2268                # update atoms in place
2269                rd.Phase['RBcoords'][:] = np.inner(trans,rd.Phase['RBcoords']).T
2270                grid.ForceRefresh()
2271                UpdateDraw()
2272
2273            def onAddVector(event):
2274                '''Adds selected atoms as a new vector rigid body.
2275                Closes out the importer tab when done.
2276                '''
2277                grid.completeEdits()
2278                rb = MakeVectorBody(os.path.split(filename)[1])
2279                UpdateVectorBody(rb,True)
2280                if len(rb['rbTypes']) < 3: return # must have at least 3 atoms
2281                rbid = ran.randint(0,sys.maxsize)
2282                data['Vector'][rbid] = rb
2283                data['RBIds']['Vector'].append(rbid)
2284                AtInfo = data['Vector']['AtInfo']
2285                for t in rb['rbTypes']:
2286                    if t in data['Vector']['AtInfo']: continue
2287                    Info = G2elem.GetAtomInfo(t)
2288                    data['Vector']['AtInfo'][t] = [Info['Drad'],Info['Color']]
2289                G2frame.G2plotNB.Delete('Rigid body')
2290                onCancel(event,0)
2291               
2292            def onAddResidue(event):
2293                '''Adds selected atoms as a new residue rigid body.
2294                Closes out the importer tab when done.
2295                '''
2296                grid.completeEdits()
2297                name = os.path.split(filename)[1]
2298                rbXYZ = []
2299                rbTypes = []
2300                atNames = []
2301                for i in rd.Phase['RBindex']:
2302                    if rd.Phase['RBselection'][i]:
2303                        rbXYZ.append(rd.Phase['RBcoords'][i])
2304                        rbTypes.append(rd.Phase['RBtypes'][i])
2305                        atNames.append(rd.Phase['RBlbls'][i])
2306                if len(rbTypes) < 3: return # must have at least 3 atoms
2307                rbXYZ = np.array(rbXYZ)
2308                rbid = ran.randint(0,sys.maxsize)
2309                data['Residue'][rbid] = {'RBname':name,'rbXYZ':rbXYZ,
2310                    'rbTypes':rbTypes,'atNames':atNames,'rbRef':[0,1,2,False],
2311                    'rbSeq':[],'SelSeq':[0,0],'useCount':0}
2312                data['RBIds']['Residue'].append(rbid)
2313                AtInfo = data['Residue']['AtInfo']
2314                for t in rbTypes:
2315                    if t in data['Residue']['AtInfo']: continue
2316                    Info = G2elem.GetAtomInfo(t)
2317                    data['Residue']['AtInfo'][t] = [Info['Drad'],Info['Color']]
2318
2319                print ('Rigid body added')
2320                G2frame.G2plotNB.Delete('Rigid body')
2321                onCancel(event,1)
2322
2323            if G2frame.rbBook.FindPage(pagename) is not None:
2324                G2frame.rbBook.DeletePage(G2frame.rbBook.FindPage(pagename))
2325            RBImp = wx.ScrolledWindow(G2frame.rbBook)
2326            RBImpPnl = wx.Panel(RBImp)
2327            G2frame.rbBook.AddPage(RBImp,pagename)
2328            G2frame.rbBook.SetSelection(G2frame.rbBook.FindPage(pagename))
2329            generalData = rd.Phase['General']
2330            AtInfo = {}
2331            ct = 1
2332            for t in rd.Phase['RBtypes']:
2333                if t in AtInfo: continue
2334                Info = G2elem.GetAtomInfo(t)
2335                AtInfo[t] = [Info['Drad'],Info['Color']]
2336            plotDefaults = {'oldxy':[0.,0.],'Quaternion':[0.,0.,0.,1.],'cameraPos':30.,'viewDir':[0,0,1],}
2337
2338            rd.Phase['RBindex'] = list(range(len(rd.Phase['RBtypes'])))
2339            rd.Phase['RBselection'] = len(rd.Phase['RBtypes']) * [1]
2340            rbData = MakeVectorBody()
2341            DrawCallback = G2plt.PlotRigidBody(G2frame,'Vector',
2342                                    AtInfo,rbData,plotDefaults)
2343
2344            mainSizer = wx.BoxSizer(wx.HORIZONTAL)
2345            btnSizer = wx.BoxSizer(wx.VERTICAL)
2346            btnSizer.Add(
2347                wx.StaticText(RBImpPnl,wx.ID_ANY,'Reorder atoms by dragging'),
2348                0,wx.ALL)
2349            btnSizer.Add((-1,15))
2350            btn = wx.Button(RBImpPnl, wx.ID_ANY, 'Set All')
2351            btn.Bind(wx.EVT_BUTTON,onSetAll)
2352            btnSizer.Add(btn,0,wx.ALIGN_CENTER)
2353            btn = wx.Button(RBImpPnl, wx.ID_ANY, 'Toggle')
2354            btn.Bind(wx.EVT_BUTTON,onToggle)
2355            btnSizer.Add(btn,0,wx.ALIGN_CENTER)
2356            btnSizer.Add((-1,15))
2357            btnSizer.Add(
2358                wx.StaticText(RBImpPnl,wx.ID_ANY,'Reorient using selected\natoms...'),
2359                0,wx.ALL)
2360            btnSizer.Add((-1,5))
2361            btn = wx.Button(RBImpPnl, wx.ID_ANY, 'Set origin')
2362            btn.Bind(wx.EVT_BUTTON,onSetOrigin)
2363            btnSizer.Add(btn,0,wx.ALIGN_CENTER)
2364
2365            bntOpts = {'plane':'xy','direction':'x'}
2366            inSizer = wx.BoxSizer(wx.HORIZONTAL)
2367            btn = wx.Button(RBImpPnl, wx.ID_ANY, 'Place in plane')
2368            btn.Bind(wx.EVT_BUTTON,onSetPlane)
2369            inSizer.Add(btn)
2370            inSizer.Add(G2G.G2ChoiceButton(RBImpPnl,('xy','yz','xz'),
2371                                           None,None,bntOpts,'plane'))
2372            btnSizer.Add(inSizer,0,wx.ALIGN_CENTER)
2373           
2374            inSizer = wx.BoxSizer(wx.HORIZONTAL)
2375            btn = wx.Button(RBImpPnl, wx.ID_ANY, 'Define as')
2376            btn.Bind(wx.EVT_BUTTON,onSetX)
2377            inSizer.Add(btn)
2378            inSizer.Add(G2G.G2ChoiceButton(RBImpPnl,('x','y','z'),
2379                                           None,None,bntOpts,'direction'))
2380            btnSizer.Add(inSizer,0,wx.ALIGN_CENTER)
2381           
2382            btnSizer.Add((-1,15))
2383            btnSizer.Add(
2384                wx.StaticText(RBImpPnl,wx.ID_ANY,'Use selected atoms to\ncreate...'),
2385                0,wx.ALL)
2386            btnSizer.Add((-1,5))
2387            btn = wx.Button(RBImpPnl, wx.ID_ANY, 'a Vector Body')
2388            btn.Bind(wx.EVT_BUTTON,onAddVector)
2389            btnSizer.Add(btn,0,wx.ALIGN_CENTER)
2390            btn = wx.Button(RBImpPnl, wx.ID_ANY, 'a Residue Body')
2391            btn.Bind(wx.EVT_BUTTON,onAddResidue)
2392            btnSizer.Add(btn,0,wx.ALIGN_CENTER)
2393            btn = wx.Button(RBImpPnl, wx.ID_CANCEL)
2394            btn.Bind(wx.EVT_BUTTON,onCancel)
2395            btnSizer.Add((-1,10))
2396            btnSizer.Add(btn,0,wx.ALIGN_CENTER)
2397
2398            mainSizer.Add(btnSizer)
2399            mainSizer.Add((5,5))
2400            grid = DragableRBGrid(RBImpPnl,rd.Phase,UpdateDraw)
2401            mainSizer.Add(grid)
2402            RBImpPnl.SetSizer(mainSizer,True)
2403            mainSizer.Layout()   
2404            Size = mainSizer.GetMinSize()
2405            Size[0] += 40
2406            Size[1] = max(Size[1],G2frame.GetSize()[1]-200) + 20
2407            RBImpPnl.SetSize(Size)
2408            RBImp.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
2409            RBImp.Scroll(0,0)
2410
2411        def GetCoords(atmsel):
2412            '''Create orthogonal coordinates for selected atoms.
2413            Place the origin at the center of the body
2414            '''
2415            atms = rd.Phase['Atoms']
2416            cell = rd.Phase['General']['Cell'][1:7]
2417            Amat,Bmat = G2lat.cell2AB(cell)
2418            rd.Phase['RBcoords'] = np.array([np.inner(Amat,atms[i][3:6]) for i in atmsel])
2419            rd.Phase['RBcoords'] -= rd.Phase['RBcoords'].mean(axis=0)  # origin to middle
2420            rd.Phase['RBtypes'] = [atms[i][1] for i in atmsel]
2421            rd.Phase['RBlbls'] = [atms[i][0] for i in atmsel]
2422
2423        def UpdateVectorBody(rb,useSelection=False):
2424            '''Put the atoms in order to pass for plotting or for storage as
2425            a vector rigid body.
2426
2427            :param dict rb: rigid body contents created in :func:`MakeVectorBody`
2428            :param bool useSelection: True if the rd.Phase['RBselection']
2429              values will be used to select which atoms are included in the
2430              rigid body. If False (default) they are included in rb
2431              and are used for plotting.         
2432            '''
2433            coordlist = []
2434            typeslist = []
2435            sellist = []
2436            for i in rd.Phase['RBindex']:
2437                use = True
2438                if useSelection and not rd.Phase['RBselection'][i]: use = False
2439                if use:
2440                    coordlist.append(rd.Phase['RBcoords'][i])
2441                    typeslist.append(rd.Phase['RBtypes'][i])
2442                    sellist.append(rd.Phase['RBselection'][i])
2443            coordlist = np.array(coordlist)
2444            rb['rbXYZ'] = coordlist
2445            rb['rbVect'] = [coordlist]
2446            rb['rbTypes'] = typeslist
2447            if not useSelection:
2448                rb['Selection'] = sellist
2449            elif 'Selection' in rb:
2450                del rb['Selection']
2451
2452        def MakeVectorBody(name=''):
2453            '''Make the basic vector rigid body dict (w/o coordinates) used for
2454            export and for plotting
2455            '''
2456            nTrans = 1
2457            vecMag = [1.0]
2458            vecRef = [False]
2459            rb = {'RBname':name,'VectMag':vecMag,
2460                    'rbRef':[0,1,2,False],'VectRef':vecRef,
2461                    'useCount':0}
2462            UpdateVectorBody(rb)
2463            return rb
2464
2465        # too lazy to figure out why wx crashes
2466        if wx.__version__.split('.')[0] != '4':
2467            wx.MessageBox('Sorry, wxPython 4.x is required to run this command',
2468                                  caption='Update Python',
2469                                  style=wx.ICON_EXCLAMATION)
2470            return
2471        if platform.python_version()[:1] == '2':
2472            wx.MessageBox('Sorry, Python >=3.x is required to run this command',
2473                                  caption='Update Python',
2474                                  style=wx.ICON_EXCLAMATION)
2475            return
2476
2477        # get importer type and a phase file of that type
2478        G2sc.LoadG2fil()
2479        choices = [rd.formatName for  rd in G2sc.Readers['Phase']] 
2480        dlg = G2G.G2SingleChoiceDialog(G2frame,'Select the format of the file',
2481                                     'select format',choices)
2482        dlg.CenterOnParent()
2483        try:
2484            if dlg.ShowModal() == wx.ID_OK:
2485                col = dlg.GetSelection()
2486            else:
2487                col = None
2488        finally:
2489            dlg.Destroy()
2490        reader = G2sc.Readers['Phase'][col]
2491
2492        choices = reader.formatName + " file ("
2493        w = ""
2494        for extn in reader.extensionlist:
2495            if w != "": w += ";"
2496            w += "*" + extn
2497        choices += w + ")|" + w
2498        #choices += "|zip archive (.zip)|*.zip"
2499        if not reader.strictExtension:
2500            choices += "|any file (*.*)|*.*"
2501        typ = '( type '+reader.formatName+')'
2502        filelist = G2G.GetImportFile(G2frame,
2503                        message="Choose phase input file"+typ,
2504                        defaultFile="",wildcard=choices,style=wx.FD_OPEN)
2505        if len(filelist) != 1: return
2506
2507        # read in the phase file
2508        filename = filelist[0]
2509        rd = reader
2510        with open(filename, 'Ur') as fp:
2511            rd.ReInitialize()
2512            rd.errors = ""
2513            if not rd.ContentsValidator(filename):   # Report error
2514                G2fl.G2Print("Warning: File {} has a validation error".format(filename))
2515                return
2516            if len(rd.selections) > 1:
2517                print("File {} has {} phases. This is unexpected."
2518                                    .format(filename,len(rd.selections)))
2519                return
2520
2521            rd.objname = os.path.basename(filename)
2522            try:
2523                flag = rd.Reader(filename)
2524            except Exception as msg:
2525                G2fl.G2Print("Warning: read of file {} failed\n{}".format(
2526                    filename,rd.errors))
2527                if GSASIIpath.GetConfigValue('debug'):
2528                    print(msg)
2529                    import traceback
2530                    print (traceback.format_exc())
2531                    GSASIIpath.IPyBreak()
2532                return
2533
2534        pagename = 'Rigid body importer'
2535        Page1()
2536        return
2537
2538    def AddVectTrans(event):
2539        'Add a translation to an existing vector rigid body'
2540        choices = []
2541        rbIdlist = []
2542        for rbid in data['RBIds']['Vector']:
2543            if rbid != 'AtInfo':
2544                rbIdlist.append(rbid)
2545                choices.append(data['Vector'][rbid]['RBname'])
2546        if len(choices) == 0:
2547            G2G.G2MessageBox(G2frame,'No Vector Rigid Bodies found',
2548                                 'No VR Bodies')
2549            return
2550        elif len(choices) == 1:
2551            rbid = rbIdlist[0]
2552        else:
2553            dlg = G2G.G2SingleChoiceDialog(G2frame,'Select the rigid body to save',
2554                                  'select format',choices)
2555            try:
2556                if dlg.ShowModal() == wx.ID_OK:
2557                    rbid = rbIdlist[dlg.GetSelection()]
2558                else:
2559                    return
2560            finally:
2561                dlg.Destroy()
2562        data['Vector'][rbid]['VectMag'] += [1.0]
2563        data['Vector'][rbid]['VectRef'] += [False]
2564        nAtoms = len(data['Vector'][rbid]['rbXYZ'])
2565        data['Vector'][rbid]['rbVect'] += [np.zeros((nAtoms,3))]
2566        UpdateVectorRB()
2567       
2568    def SaveVectorRB(event):
2569        choices = []
2570        rbIdlist = []
2571        for rbid in data['RBIds']['Vector']:
2572            if rbid != 'AtInfo':
2573                rbIdlist.append(rbid)
2574                choices.append(data['Vector'][rbid]['RBname'])
2575        if len(choices) == 0:
2576            G2G.G2MessageBox(G2frame,'No Vector Rigid Bodies found',
2577                                 'No VR Bodies')
2578            return
2579        elif len(choices) == 1:
2580            rbid = rbIdlist[0]
2581        else:
2582            dlg = G2G.G2SingleChoiceDialog(G2frame,'Select the rigid body to save',
2583                                  'select format',choices)
2584            try:
2585                if dlg.ShowModal() == wx.ID_OK:
2586                    rbid = rbIdlist[dlg.GetSelection()]
2587                else:
2588                    return
2589            finally:
2590                dlg.Destroy()
2591             
2592        pth = G2G.GetExportPath(G2frame)
2593        dlg = wx.FileDialog(G2frame, 'Choose file to save vector rigid body',
2594            pth, '', 'VRB files (*.vecbody)|*.vecbody',
2595            wx.FD_SAVE|wx.FD_OVERWRITE_PROMPT)
2596        try:
2597            if dlg.ShowModal() == wx.ID_OK:
2598                filename = dlg.GetPath()
2599                filename = os.path.splitext(filename)[0]+'.vecbody'  # set extension
2600                fp = open(filename,'w')
2601                fp.write('Name: '+data['Vector'][rbid]['RBname']+'\n')
2602                fp.write('Trans: ')
2603                for i in data['Vector'][rbid]['VectMag']:
2604                    fp.write(str(i)+" ") 
2605                fp.write('\n')
2606                ntrans = len(data['Vector'][rbid]['VectMag'])
2607                for i,sym in enumerate(data['Vector'][rbid]['rbTypes']):
2608                    fp.write("{:3s}".format(sym))
2609                    for j in range(ntrans):
2610                        fp.write('{:8.5f}{:9.5f}{:9.5f}   '
2611                            .format(*data['Vector'][rbid]['rbVect'][j][i]))
2612                    fp.write('\n')
2613                fp.close()
2614                print ('Vector rigid body saved to: '+filename)
2615        finally:
2616            dlg.Destroy()
2617           
2618    def ReadVectorRB(event):
2619        AtInfo = data['Vector']['AtInfo']
2620        pth = G2G.GetExportPath(G2frame)
2621        dlg = wx.FileDialog(G2frame, 'Choose file to read vector rigid body',
2622            pth, '', 'VRB files (*.vecbody)|*.vecbody',
2623            wx.FD_OPEN)
2624        try:
2625            if dlg.ShowModal() == wx.ID_OK:
2626                filename = dlg.GetPath()
2627                filename = os.path.splitext(filename)[0]+'.vecbody'  # set extension
2628                fp = open(filename,'r')
2629                l = fp.readline().strip()
2630                if 'Name' not in l:
2631                    fp.close()
2632                    G2frame.ErrorDialog('Read Error',
2633                        'File '+filename+' does not start with Name\nFirst line ='
2634                        +l+'\ninvalid file',parent=G2frame)
2635                    return
2636                name = l.split(':')[1].strip()
2637                trans = fp.readline().strip().split(':')[1].split()
2638                vecMag = [float(i) for i in trans]
2639                ntrans = len(trans)
2640                vecs = [[] for i in range(ntrans)]
2641                types = []
2642                l = fp.readline().strip()
2643                while l:
2644                    nums = l.strip().split()
2645                    types.append(nums.pop(0))
2646                    t = types[-1]
2647                    if t not in AtInfo:
2648                        Info = G2elem.GetAtomInfo(t)
2649                        AtInfo[t] = [Info['Drad'],Info['Color']]
2650                    for i in range(ntrans):
2651                        vecs[i].append([float(nums.pop(0)) for j in range(3)])
2652                    l = fp.readline().strip()
2653                fp.close()
2654            else:
2655                return       
2656        finally:
2657            dlg.Destroy()
2658        natoms = len(types)
2659        vecs = [np.array(vecs[i]) for i in range(ntrans)]
2660        rbid = ran.randint(0,sys.maxsize)
2661        data['Vector'][rbid] = {'RBname':name,'VectMag':vecMag,
2662                'rbXYZ':np.zeros((natoms,3)),
2663                'rbRef':[0,1,2,False],'VectRef':ntrans*[False],
2664                'rbTypes':types,
2665                'rbVect':vecs,'useCount':0}
2666        data['RBIds']['Vector'].append(rbid)
2667        UpdateVectorRB()
2668       
2669    def AddResidueRB(event):
2670        global resRBsel
2671        AtInfo = data['Residue']['AtInfo']
2672        macro = getMacroFile('rigid body')
2673        if not macro:
2674            return
2675        macStr = macro.readline()
2676        while macStr:
2677            items = macStr.split()
2678            if 'I' == items[0]:
2679                resRBsel = ran.randint(0,sys.maxsize)
2680                rbName = items[1]
2681                rbTypes = []
2682                rbXYZ = []
2683                rbSeq = []
2684                atNames = []
2685                nAtms,nSeq,nOrig,mRef,nRef = [int(items[i]) for i in [2,3,4,5,6]]
2686                for iAtm in range(nAtms):
2687                    macStr = macro.readline().split()
2688                    atName = macStr[0]
2689                    atType = macStr[1]
2690                    atNames.append(atName)
2691                    rbXYZ.append([float(macStr[i]) for i in [2,3,4]])
2692                    rbTypes.append(atType)
2693                    if atType not in AtInfo:
2694                        Info = G2elem.GetAtomInfo(atType)
2695                        AtInfo[atType] = [Info['Drad'],Info['Color']]
2696                rbXYZ = np.array(rbXYZ)-np.array(rbXYZ[nOrig-1])
2697                for iSeq in range(nSeq):
2698                    macStr = macro.readline().split()
2699                    mSeq = int(macStr[0])
2700                    for jSeq in range(mSeq):
2701                        macStr = macro.readline().split()
2702                        iBeg = int(macStr[0])-1
2703                        iFin = int(macStr[1])-1
2704                        angle = 0.0
2705                        nMove = int(macStr[2])
2706                        iMove = [int(macStr[i])-1 for i in range(3,nMove+3)]
2707                        rbSeq.append([iBeg,iFin,angle,iMove])
2708                data['Residue'][resRBsel] = {'RBname':rbName,'rbXYZ':rbXYZ,'rbTypes':rbTypes,
2709                    'atNames':atNames,'rbRef':[nOrig-1,mRef-1,nRef-1,True],'rbSeq':rbSeq,
2710                    'SelSeq':[0,0],'useCount':0}
2711                data['RBIds']['Residue'].append(resRBsel)
2712                print ('Rigid body '+rbName+' added')
2713            macStr = macro.readline()
2714        macro.close()
2715        UpdateResidueRB()
2716       
2717    def ImportResidueRB():
2718        global resRBsel
2719        AtInfo = data['Residue']['AtInfo']
2720        text,ext = getTextFile()
2721        if not text:
2722            return
2723        resRBsel = ran.randint(0,sys.maxsize)
2724        rbTypes = []
2725        rbXYZ = []
2726        atNames = []
2727        txtStr = text.readline()
2728        if 'xyz' in ext:
2729            txtStr = text.readline()
2730            txtStr = text.readline()
2731        elif 'mol2' in ext:
2732            while 'ATOM' not in txtStr:
2733                txtStr = text.readline()
2734            txtStr = text.readline()
2735        elif 'pdb' in ext:
2736            while 'ATOM' not in txtStr[:6] and 'HETATM' not in txtStr[:6]:
2737                txtStr = text.readline()
2738        items = txtStr.split()
2739        nat = 1
2740        while len(items):
2741            if 'txt' in ext:
2742                atName = items[0]
2743                atType = items[1]
2744                rbXYZ.append([float(items[i]) for i in [2,3,4]])
2745            elif 'xyz' in ext:
2746                atType = items[0]
2747                rbXYZ.append([float(items[i]) for i in [1,2,3]])
2748                atName = '%s%d'%(atType,nat)
2749            elif 'mol2' in ext:
2750                atType = items[1]
2751                atName = items[1]+items[0]
2752                rbXYZ.append([float(items[i]) for i in [2,3,4]])
2753            elif 'pdb' in ext:
2754                atType = items[-1]
2755                if not items[2][-1].isnumeric():
2756                    atName = '%s%d'%(items[2],nat)
2757                else:
2758                    atName = '5s'%items[2]
2759                xyz = txtStr[30:55].split()                   
2760                rbXYZ.append([float(x) for x in xyz])
2761            atNames.append(atName)
2762            rbTypes.append(atType)
2763            if atType not in AtInfo:
2764                Info = G2elem.GetAtomInfo(atType)
2765                AtInfo[atType] = [Info['Drad'],Info['Color']]
2766            txtStr = text.readline()
2767            if 'mol2' in ext and 'BOND' in txtStr:
2768                break
2769            if 'pdb' in ext and ('ATOM' not in txtStr[:6] and 'HETATM' not in txtStr[:6]):
2770                break
2771            items = txtStr.split()
2772            nat += 1
2773        if len(atNames) < 3:
2774            G2G.G2MessageBox(G2frame,'Not enough atoms in rigid body; must be 3 or more')
2775        else:
2776            rbXYZ = np.array(rbXYZ)-np.array(rbXYZ[0])
2777            Xxyz = rbXYZ[1]
2778            X = Xxyz/np.sqrt(np.sum(Xxyz**2))
2779            Yxyz = rbXYZ[2]
2780            Y = Yxyz/np.sqrt(np.sum(Yxyz**2))
2781            Mat = G2mth.getRBTransMat(X,Y)
2782            rbXYZ = np.inner(Mat,rbXYZ).T
2783            data['Residue'][resRBsel] = {'RBname':'UNKRB','rbXYZ':rbXYZ,'rbTypes':rbTypes,
2784                'atNames':atNames,'rbRef':[0,1,2,False],'rbSeq':[],'SelSeq':[0,0],'useCount':0}
2785            data['RBIds']['Residue'].append(resRBsel)
2786            print ('Rigid body UNKRB added')
2787        text.close()
2788        UpdateResidueRB()
2789       
2790    def SaveResidueRB():
2791        global resRBsel
2792        pth = G2G.GetExportPath(G2frame)
2793        dlg = wx.FileDialog(G2frame, 'Choose PDB file for Atom XYZ', pth, '', 
2794            'PDB files (*.pdb)|*.pdb',wx.FD_SAVE|wx.FD_OVERWRITE_PROMPT)
2795        try:
2796            if dlg.ShowModal() == wx.ID_OK:
2797                filename = dlg.GetPath()
2798                filename = os.path.splitext(filename)[0]+'.pdb'  # make extension .pdb
2799                File = open(filename,'w')       
2800                rbData =  data['Residue'][resRBsel]
2801                for iat,xyz in enumerate(rbData['rbXYZ']):
2802                    File.write('ATOM %6d  %-4s%3s     1    %8.3f%8.3f%8.3f  1.00  0.00          %2s\n'%(
2803                        iat,rbData['atNames'][iat],rbData['RBname'][:3],xyz[0],xyz[1],xyz[2],rbData['rbTypes'][iat]))
2804                File.close()
2805                print ('Atom XYZ saved to: '+filename)
2806        finally:
2807            dlg.Destroy()
2808       
2809       
2810    def FindNeighbors(Orig,XYZ,atTypes,atNames,AtInfo):
2811        Radii = []
2812        for Atype in atTypes:
2813            Radii.append(AtInfo[Atype][0])
2814        Radii = np.array(Radii)
2815        Neigh = []
2816        Dx = XYZ-XYZ[Orig]
2817        dist = np.sqrt(np.sum(Dx**2,axis=1))
2818        sumR = Radii[Orig]+Radii
2819        IndB = ma.nonzero(ma.masked_greater(dist-0.85*sumR,0.))
2820        for j in IndB[0]:
2821            if j != Orig and atTypes[j] != 'H':
2822                Neigh.append(atNames[j])
2823        return Neigh
2824       
2825    def FindAllNeighbors(XYZ,atTypes,atNames,AtInfo):
2826        NeighDict = {}
2827        for iat,xyz in enumerate(atNames):
2828            NeighDict[atNames[iat]] = FindNeighbors(iat,XYZ,atTypes,atNames,AtInfo)
2829        return NeighDict
2830       
2831    def FindRiding(Orig,Pivot,NeighDict):
2832        riding = [Orig,Pivot]
2833        iAdd = 1
2834        new = True
2835        while new:
2836            newAtms = NeighDict[riding[iAdd]]
2837            for At in newAtms:
2838                new = False
2839                if At not in riding:
2840                    riding.append(At)
2841                    new = True
2842            iAdd += 1
2843            if iAdd < len(riding):
2844                new = True
2845        return riding[2:]
2846                       
2847    def OnDefineTorsSeq(event):
2848        global resRBsel
2849        rbData = data['Residue'][resRBsel]
2850        if not len(rbData):
2851            return
2852        atNames = rbData['atNames']
2853        AtInfo = data['Residue']['AtInfo']
2854        atTypes = rbData['rbTypes']
2855        XYZ = rbData['rbXYZ']
2856        neighDict = FindAllNeighbors(XYZ,atTypes,atNames,AtInfo)
2857        TargList = []           
2858        dlg = wx.SingleChoiceDialog(G2frame,'Select origin atom for torsion sequence','Origin atom',rbData['atNames'])
2859        if dlg.ShowModal() == wx.ID_OK:
2860            Orig = dlg.GetSelection()
2861            TargList = neighDict[atNames[Orig]]
2862        dlg.Destroy()
2863        if not len(TargList):
2864            return
2865        dlg = wx.SingleChoiceDialog(G2frame,'Select pivot atom for torsion sequence','Pivot atom',TargList)
2866        if dlg.ShowModal() == wx.ID_OK:
2867            Piv = atNames.index(TargList[dlg.GetSelection()])
2868            riding = FindRiding(atNames[Orig],atNames[Piv],neighDict)
2869            Riding = []
2870            for atm in riding:
2871                Riding.append(atNames.index(atm))
2872            rbData['rbSeq'].append([Orig,Piv,0.0,Riding])           
2873        dlg.Destroy()
2874        UpdateResidueRB()
2875
2876    def UpdateVectorRB(Scroll=0):
2877        AtInfo = data['Vector']['AtInfo']
2878        refChoice = {}
2879        if 'DELETED' in str(G2frame.GetStatusBar()):   #seems to be no other way to do this (wx bug)
2880            if GSASIIpath.GetConfigValue('debug'):
2881                print ('DBG_wx error: Rigid Body/Status not cleanly deleted after Refine')
2882            return
2883        SetStatusLine(' You may use e.g. "c60" or "s60" for a vector entry')
2884        def rbNameSizer(rbid,rbData):
2885
2886            def OnRBName(event):
2887                event.Skip()
2888                Obj = event.GetEventObject()
2889                rbData['RBname'] = Obj.GetValue()
2890               
2891            def OnDelRB(event):
2892                Obj = event.GetEventObject()
2893                rbid = Indx[Obj.GetId()]
2894                if rbid in data['Vector']:
2895                    del data['Vector'][rbid]
2896                    data['RBIds']['Vector'].remove(rbid)
2897                    rbData['useCount'] -= 1
2898                wx.CallAfter(UpdateVectorRB)
2899            def OnPlotRB(event):
2900                Obj = event.GetEventObject()
2901                G2plt.PlotRigidBody(G2frame,'Vector',AtInfo,rbData,plotDefaults)
2902           
2903            nameSizer = wx.BoxSizer(wx.HORIZONTAL)
2904            nameSizer.Add(wx.StaticText(VectorRBDisplay,-1,'Rigid body name: '),
2905                0,wx.ALIGN_CENTER_VERTICAL)
2906            RBname = wx.TextCtrl(VectorRBDisplay,-1,rbData['RBname'])
2907            Indx[RBname.GetId()] = rbid
2908            RBname.Bind(wx.EVT_TEXT_ENTER,OnRBName)
2909            RBname.Bind(wx.EVT_KILL_FOCUS,OnRBName)
2910            nameSizer.Add(RBname,0,wx.ALIGN_CENTER_VERTICAL)
2911            nameSizer.Add((5,0),)
2912            plotRB =  wx.Button(VectorRBDisplay,wx.ID_ANY,'Plot',
2913                                style=wx.BU_EXACTFIT)
2914            plotRB.Bind(wx.EVT_BUTTON, OnPlotRB)
2915            Indx[plotRB.GetId()] = rbid
2916            nameSizer.Add(plotRB,0,wx.ALIGN_CENTER_VERTICAL)
2917            nameSizer.Add((5,0),)
2918            if not rbData['useCount']:
2919                delRB = wx.Button(VectorRBDisplay,wx.ID_ANY,"Delete",
2920                                style=wx.BU_EXACTFIT)
2921                delRB.Bind(wx.EVT_BUTTON, OnDelRB)
2922                Indx[delRB.GetId()] = rbid
2923                nameSizer.Add(delRB,0,wx.ALIGN_CENTER_VERTICAL)
2924            return nameSizer
2925           
2926        def rbRefAtmSizer(rbid,rbData):
2927           
2928            def OnRefSel(event):
2929                Obj = event.GetEventObject()
2930                iref = Indx[Obj.GetId()]
2931                sel = Obj.GetValue()
2932                rbData['rbRef'][iref] = atNames.index(sel)
2933                FillRefChoice(rbid,rbData)
2934           
2935            refAtmSizer = wx.BoxSizer(wx.HORIZONTAL)
2936            atNames = [name+str(i) for i,name in enumerate(rbData['rbTypes'])]
2937            rbRef = rbData.get('rbRef',[0,1,2,False])
2938            rbData['rbRef'] = rbRef
2939            if rbData['useCount']:
2940                refAtmSizer.Add(wx.StaticText(VectorRBDisplay,-1,
2941                    'Orientation reference atoms A-B-C: %s, %s, %s'%(atNames[rbRef[0]], \
2942                     atNames[rbRef[1]],atNames[rbRef[2]])),0)
2943            else:
2944                refAtmSizer.Add(wx.StaticText(VectorRBDisplay,-1,
2945                    'Orientation reference atoms A-B-C: '),0,wx.ALIGN_CENTER_VERTICAL)
2946                for i in range(3):
2947                    choices = [atNames[j] for j in refChoice[rbid][i]]
2948                    refSel = wx.ComboBox(VectorRBDisplay,-1,value='',
2949                        choices=choices,style=wx.CB_READONLY|wx.CB_DROPDOWN)
2950                    refSel.SetValue(atNames[rbRef[i]])
2951                    refSel.Bind(wx.EVT_COMBOBOX, OnRefSel)
2952                    Indx[refSel.GetId()] = i
2953                    refAtmSizer.Add(refSel,0,wx.ALIGN_CENTER_VERTICAL)
2954            return refAtmSizer
2955                       
2956        def rbVectMag(rbid,imag,rbData):
2957           
2958            def OnRBVectorMag(event):
2959                event.Skip()
2960                Obj = event.GetEventObject()
2961                rbid,imag = Indx[Obj.GetId()]
2962                try:
2963                    val = float(Obj.GetValue())
2964                    if val <= 0.:
2965                        raise ValueError
2966                    rbData['VectMag'][imag] = val
2967                except ValueError:
2968                    pass
2969                Obj.SetValue('%8.4f'%(val))
2970                wx.CallAfter(UpdateVectorRB,VectorRB.GetScrollPos(wx.VERTICAL))
2971                G2plt.PlotRigidBody(G2frame,'Vector',AtInfo,data['Vector'][rbid],plotDefaults)
2972               
2973            def OnRBVectorRef(event):
2974                Obj = event.GetEventObject()
2975                rbid,imag = Indx[Obj.GetId()]
2976                rbData['VectRef'][imag] = Obj.GetValue()
2977                       
2978            magSizer = wx.BoxSizer(wx.HORIZONTAL)
2979            magSizer.Add(wx.StaticText(VectorRBDisplay,-1,'Translation magnitude: '),
2980                0,wx.ALIGN_CENTER_VERTICAL)
2981            magValue = wx.TextCtrl(VectorRBDisplay,-1,'%8.4f'%(rbData['VectMag'][imag]))
2982            Indx[magValue.GetId()] = [rbid,imag]
2983            magValue.Bind(wx.EVT_TEXT_ENTER,OnRBVectorMag)
2984            magValue.Bind(wx.EVT_KILL_FOCUS,OnRBVectorMag)
2985            magSizer.Add(magValue,0,wx.ALIGN_CENTER_VERTICAL)
2986            magSizer.Add((5,0),)
2987            magref = wx.CheckBox(VectorRBDisplay,-1,label=' Refine?') 
2988            magref.SetValue(rbData['VectRef'][imag])
2989            magref.Bind(wx.EVT_CHECKBOX,OnRBVectorRef)
2990            Indx[magref.GetId()] = [rbid,imag]
2991            magSizer.Add(magref,0,wx.ALIGN_CENTER_VERTICAL)
2992            return magSizer
2993           
2994        def rbVectors(rbid,imag,mag,XYZ,rbData):
2995
2996            def TypeSelect(event):
2997                AtInfo = data['Vector']['AtInfo']
2998                r,c = event.GetRow(),event.GetCol()
2999                if vecGrid.GetColLabelValue(c) == 'Type':
3000                    PE = G2elemGUI.PickElement(G2frame,oneOnly=True)
3001                    if PE.ShowModal() == wx.ID_OK:
3002                        if PE.Elem != 'None':
3003                            El = PE.Elem.strip().lower().capitalize()
3004                            if El not in AtInfo:
3005                                Info = G2elem.GetAtomInfo(El)
3006                                AtInfo[El] = [Info['Drad'],Info['Color']]
3007                            rbData['rbTypes'][r] = El
3008                            vecGrid.SetCellValue(r,c,El)
3009                    PE.Destroy()
3010                wx.CallAfter(UpdateVectorRB,VectorRB.GetScrollPos(wx.VERTICAL))
3011
3012            def ChangeCell(event):
3013                r,c =  event.GetRow(),event.GetCol()
3014                if r >= 0 and (0 <= c < 3):
3015                    try:
3016                        val = float(vecGrid.GetCellValue(r,c))
3017                        rbData['rbVect'][imag][r][c] = val
3018                    except ValueError:
3019                        pass
3020                G2plt.PlotRigidBody(G2frame,'Vector',AtInfo,data['Vector'][rbid],plotDefaults)
3021                wx.CallAfter(UpdateVectorRB,VectorRB.GetScrollPos(wx.VERTICAL))
3022
3023            vecSizer = wx.BoxSizer()
3024            Types = 3*[wg.GRID_VALUE_FLOAT+':10,5',]+[wg.GRID_VALUE_STRING,]+3*[wg.GRID_VALUE_FLOAT+':10,5',]
3025            colLabels = ['Vector x','Vector y','Vector z','Type','Cart x','Cart y','Cart z']
3026            table = []
3027            rowLabels = []
3028            for ivec,xyz in enumerate(rbData['rbVect'][imag]):
3029                table.append(list(xyz)+[rbData['rbTypes'][ivec],]+list(XYZ[ivec]))
3030                rowLabels.append(str(ivec))
3031            vecTable = G2G.Table(table,rowLabels=rowLabels,colLabels=colLabels,types=Types)
3032            vecGrid = G2G.GSGrid(VectorRBDisplay)
3033            vecGrid.SetTable(vecTable, True)
3034            if 'phoenix' in wx.version():
3035                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGED, ChangeCell)
3036            else:
3037                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGE, ChangeCell)
3038            if not imag:
3039                vecGrid.Bind(wg.EVT_GRID_CELL_LEFT_DCLICK, TypeSelect)
3040            attr = wx.grid.GridCellAttr()
3041            attr.IncRef()
3042            attr.SetEditor(G2G.GridFractionEditor(vecGrid))
3043            for c in range(3):
3044                attr.IncRef()
3045                vecGrid.SetColAttr(c, attr)
3046            for row in range(vecTable.GetNumberRows()):
3047                if imag:
3048                    vecGrid.SetCellStyle(row,3,VERY_LIGHT_GREY,True)                   
3049                for col in [4,5,6]:
3050                    vecGrid.SetCellStyle(row,col,VERY_LIGHT_GREY,True)
3051#            vecGrid.SetScrollRate(0,0)
3052            vecGrid.AutoSizeColumns(False)
3053            vecSizer.Add(vecGrid)
3054            return vecSizer
3055       
3056        def FillRefChoice(rbid,rbData):
3057            choiceIds = [i for i in range(len(rbData['rbTypes']))]
3058           
3059            rbRef = rbData.get('rbRef',[-1,-1,-1,False])
3060            for i in range(3):
3061                choiceIds.remove(rbRef[i])
3062            refChoice[rbid] = [choiceIds[:],choiceIds[:],choiceIds[:]]
3063            for i in range(3):
3064                refChoice[rbid][i].append(rbRef[i])
3065                refChoice[rbid][i].sort()     
3066           
3067        if VectorRB.GetSizer(): VectorRB.GetSizer().Clear(True)
3068        VectorRBSizer = wx.BoxSizer(wx.VERTICAL)
3069        first = True
3070        for rbid in data['RBIds']['Vector']:
3071            if rbid != 'AtInfo':
3072                rbData = data['Vector'][rbid]
3073                FillRefChoice(rbid,rbData)
3074                if not first:
3075                    G2G.HorizontalLine(VectorRBSizer,VectorRBDisplay)
3076                VectorRBSizer.Add(rbNameSizer(rbid,rbData),0)
3077                VectorRBSizer.Add(rbRefAtmSizer(rbid,rbData),0)
3078                XYZ = np.array([[0.,0.,0.] for Ty in rbData['rbTypes']])
3079                for imag,mag in enumerate(rbData['VectMag']):
3080                    XYZ += mag*rbData['rbVect'][imag]
3081                    VectorRBSizer.Add(rbVectMag(rbid,imag,rbData),0)
3082                    VectorRBSizer.Add(rbVectors(rbid,imag,mag,XYZ,rbData),0)
3083                VectorRBSizer.Add((5,5),0)
3084                data['Vector'][rbid]['rbXYZ'] = XYZ       
3085                first = False
3086        VectorRBSizer.Layout()   
3087        VectorRBDisplay.SetSizer(VectorRBSizer,True)
3088        Size = VectorRBSizer.GetMinSize()
3089        Size[0] += 40
3090        Size[1] = max(Size[1],450) + 20
3091        VectorRBDisplay.SetSize(Size)
3092        VectorRB.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
3093        VectorRB.Scroll(0,Scroll)
3094       
3095    def UpdateResidueRB():
3096        '''Draw the contents of the Residue Rigid Body tab for Rigid Bodies tree entry
3097        '''
3098        global resRBsel
3099        def rbNameSizer(rbid,rbData):
3100
3101            def OnDelRB(event):
3102                Obj = event.GetEventObject()
3103                rbid = Indx[Obj.GetId()]
3104                if rbid in data['Residue']: 
3105                    del data['Residue'][rbid]
3106                    data['RBIds']['Residue'].remove(rbid)
3107                wx.CallAfter(UpdateResidueRB)
3108               
3109            def OnStripH(event):
3110                Obj = event.GetEventObject()
3111                rbid = Indx[Obj.GetId()]
3112                if rbid in data['Residue']:
3113                    newNames = []
3114                    newTypes = []
3115                    newXYZ = []
3116                    for i,atype in enumerate(rbData['rbTypes']):
3117                        if atype != 'H':
3118                            newNames.append(rbData['atNames'][i])
3119                            newTypes.append(rbData['rbTypes'][i])
3120                            newXYZ.append(rbData['rbXYZ'][i])
3121                    rbData['atNames'] = newNames
3122                    rbData['rbTypes'] = newTypes
3123                    rbData['rbXYZ'] = newXYZ
3124                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
3125                wx.CallAfter(UpdateResidueRB)
3126                   
3127            def OnPlotRB(event):
3128                Obj = event.GetEventObject()
3129                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
3130               
3131            # start of rbNameSizer
3132            nameSizer = wx.BoxSizer(wx.HORIZONTAL)
3133            nameSizer.Add(wx.StaticText(ResidueRBDisplay,-1,'Residue name: '),
3134                0,wx.ALIGN_CENTER_VERTICAL)
3135            nameSizer.Add(G2G.ValidatedTxtCtrl(ResidueRBDisplay,rbData,'RBname'),0,WACV)
3136            nameSizer.Add((5,0),)
3137            plotRB =  wx.Button(ResidueRBDisplay,wx.ID_ANY,'Plot',
3138                                style=wx.BU_EXACTFIT)
3139            plotRB.Bind(wx.EVT_BUTTON, OnPlotRB)
3140            Indx[plotRB.GetId()] = rbid
3141            nameSizer.Add(plotRB,0,wx.ALIGN_CENTER_VERTICAL)
3142            nameSizer.Add((5,0),)
3143            if not rbData['useCount']:
3144                delRB = wx.Button(ResidueRBDisplay,wx.ID_ANY,"Delete",
3145                                style=wx.BU_EXACTFIT)
3146                delRB.Bind(wx.EVT_BUTTON, OnDelRB)
3147                Indx[delRB.GetId()] = rbid
3148                nameSizer.Add(delRB,0,wx.ALIGN_CENTER_VERTICAL)
3149                if 'H'  in rbData['rbTypes']:
3150                    stripH = wx.Button(ResidueRBDisplay,wx.ID_ANY,
3151                                "Strip H-atoms",
3152                                style=wx.BU_EXACTFIT)
3153                    stripH.Bind(wx.EVT_BUTTON, OnStripH)
3154                    Indx[stripH.GetId()] = rbid
3155                    nameSizer.Add(stripH,0,wx.ALIGN_CENTER_VERTICAL)
3156            return nameSizer
3157           
3158        def rbResidues(rbid,rbData):
3159           
3160            def TypeSelect(event):
3161                AtInfo = data['Residue']['AtInfo']
3162                r,c = event.GetRow(),event.GetCol()
3163                if resGrid.GetColLabelValue(c) == 'Type':
3164                    PE = G2elemGUI.PickElement(G2frame,oneOnly=True)
3165                    if PE.ShowModal() == wx.ID_OK:
3166                        if PE.Elem != 'None':
3167                            El = PE.Elem.strip().lower().capitalize()
3168                            if El not in AtInfo:
3169                                Info = G2elem.GetAtomInfo(El)
3170                                AtInfo[El] = [Info['Drad']['Color']]
3171                            rbData['rbTypes'][r] = El
3172                            resGrid.SetCellValue(r,c,El)
3173                    PE.Destroy()
3174
3175            def ChangeCell(event):
3176                r,c =  event.GetRow(),event.GetCol()
3177                if r >= 0 and (0 <= c < 3):
3178                    try:
3179                        val = float(resGrid.GetCellValue(r,c))
3180                        rbData['rbXYZ'][r][c] = val
3181                    except ValueError:
3182                        pass
3183                       
3184            def RowSelect(event):
3185                r,c =  event.GetRow(),event.GetCol()
3186                if c < 0:                   #only row clicks
3187                    for iGrid in resList:
3188                        iGrid.ClearSelection()
3189                    resGrid.SelectRow(r,True)
3190
3191            def OnRefSel(event):
3192                Obj = event.GetEventObject()
3193                iref,res,jref = Indx[Obj.GetId()]
3194                sel = Obj.GetValue()
3195                ind = atNames.index(sel)
3196                if rbData['rbTypes'][ind] == 'H':
3197                    G2G.G2MessageBox(G2frame,'You should not select an H-atom for rigid body orientation')
3198                rbData['rbRef'][iref] = ind
3199                FillRefChoice(rbid,rbData)
3200                for i,ref in enumerate(RefObjs[jref]):
3201                    ref.SetItems([atNames[j] for j in refChoice[rbid][i]])
3202                    ref.SetValue(atNames[rbData['rbRef'][i]])                   
3203                rbXYZ = rbData['rbXYZ']
3204                if not iref:     #origin change
3205                    rbXYZ -= rbXYZ[ind]
3206                #TODO - transform all atom XYZ by axis choices
3207                Xxyz = rbXYZ[rbData['rbRef'][1]]
3208                X = Xxyz/np.sqrt(np.sum(Xxyz**2))
3209                Yxyz = rbXYZ[rbData['rbRef'][2]]
3210                Y = Yxyz/np.sqrt(np.sum(Yxyz**2))
3211                Mat = G2mth.getRBTransMat(X,Y)
3212                rbXYZ = np.inner(Mat,rbXYZ).T
3213                rbData['rbXYZ'] = rbXYZ
3214                res.ClearSelection()
3215                resTable = res.GetTable()
3216                for r in range(res.GetNumberRows()):
3217                    row = resTable.GetRowValues(r)
3218                    row[2:4] = rbXYZ[r]
3219                    resTable.SetRowValues(r,row)
3220                res.ForceRefresh()
3221                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
3222               
3223            Types = 2*[wg.GRID_VALUE_STRING,]+3*[wg.GRID_VALUE_FLOAT+':10,5',]
3224            colLabels = ['Name','Type','Cart x','Cart y','Cart z']
3225            table = []
3226            rowLabels = []
3227            for ivec,xyz in enumerate(rbData['rbXYZ']):
3228                table.append([rbData['atNames'][ivec],]+[rbData['rbTypes'][ivec],]+list(xyz))
3229                rowLabels.append(str(ivec))
3230            vecTable = G2G.Table(table,rowLabels=rowLabels,colLabels=colLabels,types=Types)
3231            resGrid = G2G.GSGrid(ResidueRBDisplay)
3232            Indx[resGrid.GetId()] = rbid
3233            resList.append(resGrid)
3234            resGrid.SetTable(vecTable, True)
3235            if 'phoenix' in wx.version():
3236                resGrid.Bind(wg.EVT_GRID_CELL_CHANGED, ChangeCell)
3237            else:
3238                resGrid.Bind(wg.EVT_GRID_CELL_CHANGE, ChangeCell)
3239            resGrid.Bind(wg.EVT_GRID_CELL_LEFT_DCLICK, TypeSelect)
3240            resGrid.Bind(wg.EVT_GRID_LABEL_LEFT_CLICK, RowSelect)
3241            attr = wx.grid.GridCellAttr()
3242            attr.IncRef()
3243            attr.SetEditor(G2G.GridFractionEditor(resGrid))
3244            for c in range(3):
3245                attr.IncRef()
3246                resGrid.SetColAttr(c, attr)
3247            for row in range(vecTable.GetNumberRows()):
3248                for col in range(5):
3249                    resGrid.SetCellStyle(row,col,VERY_LIGHT_GREY,True)
3250            resGrid.AutoSizeColumns(False)
3251            vecSizer = wx.BoxSizer()
3252            vecSizer.Add(resGrid)
3253           
3254            refAtmSizer = wx.BoxSizer(wx.HORIZONTAL)
3255            atNames = rbData['atNames']
3256            rbRef = rbData['rbRef']
3257            if rbData['rbRef'][3] or rbData['useCount']:
3258                refAtmSizer.Add(wx.StaticText(ResidueRBDisplay,-1,
3259                    'Orientation reference non-H atoms A-B-C: %s, %s, %s'%(atNames[rbRef[0]], \
3260                     atNames[rbRef[1]],atNames[rbRef[2]])),0)
3261            else:
3262                refAtmSizer.Add(wx.StaticText(ResidueRBDisplay,-1,
3263                    'Orientation reference non-H atoms A-B-C: '),0,wx.ALIGN_CENTER_VERTICAL)
3264                refObj = [0,0,0]
3265                for i in range(3):
3266                    choices = [atNames[j] for j in refChoice[rbid][i]]
3267                    refSel = wx.ComboBox(ResidueRBDisplay,-1,value='',
3268                        choices=choices,style=wx.CB_READONLY|wx.CB_DROPDOWN)
3269                    refSel.SetValue(atNames[rbRef[i]])
3270                    refSel.Bind(wx.EVT_COMBOBOX, OnRefSel)
3271                    Indx[refSel.GetId()] = [i,resGrid,len(RefObjs)]
3272                    refObj[i] = refSel
3273                    refAtmSizer.Add(refSel,0,wx.ALIGN_CENTER_VERTICAL)
3274                RefObjs.append(refObj)
3275           
3276            mainSizer = wx.BoxSizer(wx.VERTICAL)
3277            mainSizer.Add(refAtmSizer)
3278            mainSizer.Add(vecSizer)
3279            return mainSizer
3280           
3281        def Add2SeqSizer(seqSizer,angSlide,rbid,iSeq,Seq,atNames):
3282           
3283            def ChangeAngle(event):
3284                event.Skip()
3285                Obj = event.GetEventObject()
3286                rbid,Seq = Indx[Obj.GetId()][:2]
3287                val = Seq[2]
3288                try:
3289                    val = float(Obj.GetValue())
3290                    Seq[2] = val
3291                except ValueError:
3292                    pass
3293                Obj.SetValue('%8.2f'%(val))
3294                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,data['Residue'][rbid],plotDefaults)
3295               
3296            def OnRadBtn(event):
3297                Obj = event.GetEventObject()
3298                Seq,iSeq,angId = Indx[Obj.GetId()]
3299                data['Residue'][rbid]['SelSeq'] = [iSeq,angId]
3300                angSlide.SetValue(int(100*Seq[2]))
3301               
3302            def OnDelBtn(event):
3303                Obj = event.GetEventObject()
3304                rbid,Seq = Indx[Obj.GetId()]
3305                data['Residue'][rbid]['rbSeq'].remove(Seq)       
3306                wx.CallAfter(UpdateResidueRB)
3307           
3308            iBeg,iFin,angle,iMove = Seq
3309            ang = wx.TextCtrl(ResidueRBDisplay,wx.ID_ANY,
3310                    '%8.2f'%(angle),size=(70,-1),style=wx.TE_PROCESS_ENTER)
3311            if not iSeq:
3312                radBt = wx.RadioButton(ResidueRBDisplay,wx.ID_ANY,
3313                                           '',style=wx.RB_GROUP)
3314                data['Residue'][rbid]['SelSeq'] = [iSeq,ang.GetId()]
3315                radBt.SetValue(True)
3316            else:
3317                radBt = wx.RadioButton(ResidueRBDisplay,wx.ID_ANY,'')
3318            radBt.Bind(wx.EVT_RADIOBUTTON,OnRadBtn)                   
3319            seqSizer.Add(radBt)
3320            delBt =  wx.Button(ResidueRBDisplay,wx.ID_ANY,'Del',
3321                                style=wx.BU_EXACTFIT)
3322            delBt.Bind(wx.EVT_BUTTON,OnDelBtn)
3323            seqSizer.Add(delBt)
3324            bond = wx.StaticText(ResidueRBDisplay,wx.ID_ANY,
3325                        '%s %s'%(atNames[iBeg],atNames[iFin]),size=(50,20))
3326            seqSizer.Add(bond,0,wx.ALIGN_CENTER_VERTICAL)
3327            Indx[radBt.GetId()] = [Seq,iSeq,ang.GetId()]
3328            Indx[delBt.GetId()] = [rbid,Seq]
3329            Indx[ang.GetId()] = [rbid,Seq,ang]
3330            ang.Bind(wx.EVT_TEXT_ENTER,ChangeAngle)
3331            ang.Bind(wx.EVT_KILL_FOCUS,ChangeAngle)
3332            seqSizer.Add(ang,0,wx.ALIGN_CENTER_VERTICAL)
3333            atms = ''
3334            for i in iMove:   
3335                atms += ' %s,'%(atNames[i])
3336            moves = wx.StaticText(ResidueRBDisplay,wx.ID_ANY,
3337                            atms[:-1],size=(200,20))
3338            seqSizer.Add(moves,1,wx.ALIGN_CENTER_VERTICAL|wx.EXPAND|wx.RIGHT)
3339            return seqSizer
3340           
3341        def SlideSizer():
3342           
3343            def OnSlider(event):
3344                Obj = event.GetEventObject()
3345                rbData = Indx[Obj.GetId()]
3346                iSeq,angId = rbData['SelSeq']
3347                val = float(Obj.GetValue())/100.
3348                rbData['rbSeq'][iSeq][2] = val
3349                Indx[angId][2].SetValue('%8.2f'%(val))
3350                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
3351           
3352            slideSizer = wx.BoxSizer(wx.HORIZONTAL)
3353            slideSizer.Add(wx.StaticText(ResidueRBDisplay,-1,'Selected torsion angle:'),0)
3354            iSeq,angId = rbData['SelSeq']
3355            angSlide = wx.Slider(ResidueRBDisplay,-1,
3356                int(100*rbData['rbSeq'][iSeq][2]),0,36000,size=(200,20),
3357                style=wx.SL_HORIZONTAL)
3358            angSlide.Bind(wx.EVT_SLIDER, OnSlider)
3359            Indx[angSlide.GetId()] = rbData
3360            slideSizer.Add(angSlide,0)           
3361            return slideSizer,angSlide
3362           
3363        def FillRefChoice(rbid,rbData):
3364            choiceIds = [i for i in range(len(rbData['atNames']))]
3365            for seq in rbData['rbSeq']:
3366                for i in seq[3]:
3367                    try:
3368                        choiceIds.remove(i)
3369                    except ValueError:
3370                        pass
3371            rbRef = rbData['rbRef']
3372            for i in range(3):
3373                try:
3374                    choiceIds.remove(rbRef[i])
3375                except ValueError:
3376                    pass
3377            refChoice[rbid] = [choiceIds[:],choiceIds[:],choiceIds[:]]
3378            for i in range(3):
3379                refChoice[rbid][i].append(rbRef[i])
3380                refChoice[rbid][i].sort()
3381               
3382        def OnSelect(event):
3383            global resRBsel
3384            sel = select.GetSelection()
3385            if sel == 0: return # 1st entry is blank
3386            rbname = rbchoice[sel-1]
3387            resRBsel = RBnames[rbname]
3388            wx.CallLater(100,UpdateResidueRB)
3389           
3390        #----- beginning of UpdateResidueRB -----------------------------------------------
3391        AtInfo = data['Residue']['AtInfo']
3392        refChoice = {}
3393        RefObjs = []
3394
3395        GS = ResidueRBDisplay.GetSizer()
3396        if GS: 
3397            try:        #get around a c++ error in wx 4.0; doing is again seems to be OK
3398                GS.Clear(True)
3399            except:
3400                GS.Clear(True)
3401       
3402        RBnames = {}
3403        for rbid in data['RBIds']['Residue']:
3404            RBnames.update({data['Residue'][rbid]['RBname']:rbid,})
3405        if not RBnames:
3406            return
3407        rbchoice = list(RBnames.keys())
3408        if GS: 
3409            ResidueRBSizer = GS
3410        else:
3411            ResidueRBSizer = wx.BoxSizer(wx.VERTICAL)
3412        if len(RBnames) > 1:
3413            selSizer = wx.BoxSizer(wx.HORIZONTAL)
3414            selSizer.Add(wx.StaticText(ResidueRBDisplay,
3415                                label=' Select residue to view:'),0)
3416            rbchoice.sort()
3417            select = wx.ComboBox(ResidueRBDisplay,choices=['']+rbchoice)
3418            select.Bind(wx.EVT_COMBOBOX,OnSelect)
3419            selSizer.Add(select,0)
3420            ResidueRBSizer.Add(selSizer,0)
3421        if resRBsel not in data['RBIds']['Residue']:
3422            resRBsel = RBnames[rbchoice[0]]
3423        rbData = data['Residue'][resRBsel]
3424        FillRefChoice(resRBsel,rbData)
3425        ResidueRBSizer.Add(rbNameSizer(resRBsel,rbData),0)
3426        ResidueRBSizer.Add(rbResidues(resRBsel,rbData),0)
3427        if len(rbData['rbSeq']):
3428            ResidueRBSizer.Add((-1,15),0)
3429            slideSizer,angSlide = SlideSizer()
3430            seqSizer = wx.FlexGridSizer(0,5,4,8)
3431            for lbl in 'Sel','','Bond','Angle','Riding atoms':
3432                seqSizer.Add(wx.StaticText(ResidueRBDisplay,wx.ID_ANY,lbl))
3433            ResidueRBSizer.Add(seqSizer)
3434#            for iSeq,Seq in enumerate(rbData['rbSeq']):
3435#                ResidueRBSizer.Add(SeqSizer(angSlide,resRBsel,iSeq,Seq,rbData['atNames']))
3436            for iSeq,Seq in enumerate(rbData['rbSeq']):
3437                Add2SeqSizer(seqSizer,angSlide,resRBsel,iSeq,Seq,rbData['atNames'])
3438            ResidueRBSizer.Add(slideSizer,)
3439
3440        ResidueRBSizer.Add((5,25),)
3441        ResidueRBSizer.Layout()   
3442        ResidueRBDisplay.SetSizer(ResidueRBSizer,True)
3443        ResidueRBDisplay.SetAutoLayout(True)
3444        Size = ResidueRBSizer.GetMinSize()
3445        ResidueRBDisplay.SetSize(Size)
3446       
3447        Size[0] += 40
3448        Size[1] = max(Size[1],450) + 20
3449        ResidueRB.SetSize(Size)
3450        ResidueRB.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
3451        G2frame.dataWindow.SendSizeEvent()
3452       
3453        ResidueRBDisplay.Show()
3454       
3455    def SetStatusLine(text):
3456        G2frame.GetStatusBar().SetStatusText(text,1)                                     
3457
3458    G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.RigidBodyMenu)
3459    SetStatusLine('')
3460    UpdateVectorRB()
3461    G2frame.rbBook.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, OnPageChanged)
3462    wx.CallAfter(OnPageChanged,None)
3463
3464def ShowIsoDistortCalc(G2frame,phase=None):
3465    '''Compute the ISODISTORT mode values from the current coordinates.
3466    Called in response to the (Phase/Atoms tab) AtomCompute or
3467    Constraints/Edit Constr. "Show ISODISTORT modes" menu item, which
3468    should be enabled only when Phase['ISODISTORT'] is defined.
3469    '''
3470    def _onClose(event):
3471        dlg.EndModal(wx.ID_CANCEL)
3472    def fmtHelp(item,fullname):
3473        helptext = "A new variable"
3474        if item[-3]:
3475            helptext += " named "+str(item[-3])
3476        helptext += " is a linear combination of the following parameters:\n"
3477        first = True
3478        for term in item[:-3]:
3479            line = ''
3480            var = str(term[1])
3481            m = term[0]
3482            if first:
3483                first = False
3484                line += ' = '
3485            else:
3486                if m >= 0:
3487                    line += ' + '
3488                else:
3489                    line += ' - '
3490                m = abs(m)
3491            line += '%.3f*%s '%(m,var)
3492            varMean = G2obj.fmtVarDescr(var)
3493            helptext += "\n" + line + " ("+ varMean + ")"
3494        helptext += '\n\nISODISTORT full name: '+str(fullname)
3495        return helptext
3496
3497    Histograms,Phases = G2frame.GetUsedHistogramsAndPhasesfromTree() # init for constraint
3498    isophases = [p for p in Phases if 'ISODISTORT' in Phases[p]]
3499   
3500    if not isophases:
3501        G2frame.ErrorDialog('no ISODISTORT phases',
3502                            'Unexpected error: no ISODISTORT phases')
3503        return
3504    if phase and phase not in isophases:
3505        G2frame.ErrorDialog('not ISODISTORT phase',
3506                            'Unexpected error: selected phase is not an ISODISTORT phase')
3507        print('Unexpected error: selected phase is not an ISODISTORT phase',
3508                  phase,isophases)
3509    elif not phase and len(isophases) == 1:
3510        phase = isophases[0]
3511    elif not phase:
3512        dlg = wx.SingleChoiceDialog(G2frame,'Select phase from ISODISTORT phases',
3513                                        'Select Phase',isophases)
3514        if dlg.ShowModal() == wx.ID_OK:
3515            sel = dlg.GetSelection()
3516            phase = isophases[sel]
3517        else:
3518            return
3519    # if len(data.get('Histograms',[])) == 0:
3520    #     G2frame.ErrorDialog(
3521    #         'No data',
3522    #         'Sorry, this computation requires that a histogram first be added to the phase'
3523    #         )
3524    #     return
3525   
3526    covdata = G2frame.GPXtree.GetItemPyData(
3527        G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Covariance'))
3528    # make a lookup table for named NewVar Phase constraints
3529    sub = G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Constraints') 
3530    Constraints = G2frame.GPXtree.GetItemPyData(sub)
3531    constDict = {}
3532    for c in Constraints['Phase']:
3533        if c[-1] != 'f' or not c[-3]: continue
3534        constDict[str(c[-3])] = c
3535
3536    parmDict,varyList = G2frame.MakeLSParmDict()
3537
3538    dlg = wx.Dialog(G2frame,wx.ID_ANY,'ISODISTORT mode values',#size=(630,400),
3539                       style=wx.DEFAULT_DIALOG_STYLE|wx.RESIZE_BORDER)
3540    mainSizer = wx.BoxSizer(wx.VERTICAL)
3541    if 'ISODISTORT' not in Phases[phase]:
3542        G2frame.ErrorDialog('not ISODISTORT phase',
3543                            'Unexpected error: selected phase is not an ISODISTORT phase')
3544        return
3545    data = Phases[phase]
3546    ISO = data['ISODISTORT']
3547    mainSizer.Add(wx.StaticText(dlg,wx.ID_ANY,
3548                                'ISODISTORT mode computation for cordinates in phase '+
3549                                str(data['General'].get('Name'))))
3550    aSizer = wx.BoxSizer(wx.HORIZONTAL)
3551    panel1 = wxscroll.ScrolledPanel(
3552        dlg, wx.ID_ANY,#size=(100,200),
3553        style = wx.TAB_TRAVERSAL|wx.SUNKEN_BORDER)
3554    subSizer1 = wx.FlexGridSizer(cols=2,hgap=5,vgap=2)
3555    panel2 = wxscroll.ScrolledPanel(
3556        dlg, wx.ID_ANY,#size=(100,200),
3557        style = wx.TAB_TRAVERSAL|wx.SUNKEN_BORDER)
3558    subSizer2 = wx.FlexGridSizer(cols=3,hgap=5,vgap=2)
3559    subSizer1.Add(wx.StaticText(panel1,wx.ID_ANY,'Parameter name  '))
3560    subSizer1.Add(wx.StaticText(panel1,wx.ID_ANY,' value'),0,wx.ALIGN_RIGHT)
3561    subSizer2.Add((-1,-1))
3562    subSizer2.Add(wx.StaticText(panel2,wx.ID_ANY,'Mode name  '))
3563    subSizer2.Add(wx.StaticText(panel2,wx.ID_ANY,' value'),0,wx.ALIGN_RIGHT)
3564    # ISODISTORT displacive modes
3565    if 'G2VarList' in ISO:
3566        deltaList = []
3567        for gv,Ilbl in zip(ISO['G2VarList'],ISO['IsoVarList']):
3568            dvar = gv.varname()
3569            var = dvar.replace('::dA','::A')
3570            albl = Ilbl[:Ilbl.rfind('_')]
3571            v = Ilbl[Ilbl.rfind('_')+1:]
3572            pval = ISO['ParentStructure'][albl][['dx','dy','dz'].index(v)]
3573            if var in parmDict:
3574                cval = parmDict[var][0]
3575            else:
3576                dlg.EndModal(wx.ID_CANCEL)
3577                G2frame.ErrorDialog('Atom not found',"No value found for parameter "+str(var))
3578                return
3579            deltaList.append(cval-pval)
3580        modeVals = np.inner(ISO['Var2ModeMatrix'],deltaList)
3581        for lbl,xyz,var,val,norm,G2mode in zip(
3582                ISO['IsoVarList'],deltaList,
3583                ISO['IsoModeList'],modeVals,ISO['NormList'],ISO['G2ModeList'] ):
3584            #GSASIIpath.IPyBreak()
3585            if str(G2mode) in constDict:
3586                ch = G2G.HelpButton(panel2,fmtHelp(constDict[str(G2mode)],var))
3587                subSizer2.Add(ch,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)
3588            else:
3589                subSizer2.Add((-1,-1))
3590            subSizer1.Add(wx.StaticText(panel1,wx.ID_ANY,str(lbl)))
3591            try:
3592                value = G2py3.FormatSigFigs(xyz)
3593            except TypeError:
3594                value = str(xyz)           
3595            subSizer1.Add(wx.StaticText(panel1,wx.ID_ANY,value),0,wx.ALIGN_RIGHT)
3596            subSizer2.Add(wx.StaticText(panel2,wx.ID_ANY,str(var)))
3597            try:
3598                value = G2py3.FormatSigFigs(val/norm)
3599                if 'varyList' in covdata:
3600                    if str(G2mode) in covdata['varyList']:
3601                        sig = covdata['sig'][covdata['varyList'].index(str(G2mode))]
3602                        value = G2mth.ValEsd(val/norm,sig/norm)
3603            except TypeError:
3604                value = '?'
3605            subSizer2.Add(wx.StaticText(panel2,wx.ID_ANY,value),0,wx.ALIGN_RIGHT)
3606            #GSASIIpath.IPyBreak()
3607    # ISODISTORT occupancy modes
3608    if 'G2OccVarList' in ISO:
3609        deltaList = []
3610        for gv,Ilbl in zip(ISO['G2OccVarList'],ISO['OccVarList']):
3611            var = gv.varname()
3612            albl = Ilbl[:Ilbl.rfind('_')]
3613            pval = ISO['BaseOcc'][albl]
3614            if var in parmDict:
3615                cval = parmDict[var][0]
3616            else:
3617                dlg.EndModal(wx.ID_CANCEL)
3618                G2frame.ErrorDialog('Atom not found',"No value found for parameter "+str(var))
3619                return
3620            deltaList.append(cval-pval)
3621        modeVals = np.inner(ISO['Var2OccMatrix'],deltaList)
3622        for lbl,delocc,var,val,norm,G2mode in zip(
3623                ISO['OccVarList'],deltaList,
3624                ISO['OccModeList'],modeVals,ISO['OccNormList'],ISO['G2OccModeList']):
3625            if str(G2mode) in constDict:
3626                ch = G2G.HelpButton(panel2,fmtHelp(constDict[str(G2mode)],var))
3627                subSizer2.Add(ch,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)
3628            else:
3629                subSizer2.Add((-1,-1))
3630            subSizer1.Add(wx.StaticText(panel1,wx.ID_ANY,str(lbl)))
3631            try:
3632                value = G2py3.FormatSigFigs(delocc)
3633            except TypeError:
3634                value = str(delocc)
3635            subSizer1.Add(wx.StaticText(panel1,wx.ID_ANY,value),0,wx.ALIGN_RIGHT)
3636            #subSizer.Add((10,-1))
3637            subSizer2.Add(wx.StaticText(panel2,wx.ID_ANY,str(var)))
3638            try:
3639                value = G2py3.FormatSigFigs(val/norm)
3640                if 'varyList' in covdata:
3641                    if str(G2mode) in covdata['varyList']:
3642                        sig = covdata['sig'][covdata['varyList'].index(str(G2mode))]
3643                        value = G2mth.ValEsd(val/norm,sig/norm)
3644            except TypeError:
3645                value = '?'
3646            subSizer2.Add(wx.StaticText(panel2,wx.ID_ANY,value),0,wx.ALIGN_RIGHT)
3647
3648    # finish up ScrolledPanel
3649    panel1.SetSizer(subSizer1)
3650    panel2.SetSizer(subSizer2)
3651    panel1.SetAutoLayout(1)
3652    panel1.SetupScrolling()
3653    panel2.SetAutoLayout(1)
3654    panel2.SetupScrolling()
3655    # Allow window to be enlarged but not made smaller
3656    dlg.SetSizer(mainSizer)
3657    w1,l1 = subSizer1.GetSize()
3658    w2,l2 = subSizer2.GetSize()
3659    panel1.SetMinSize((w1+10,200))
3660    panel2.SetMinSize((w2+20,200))
3661    aSizer.Add(panel1,1, wx.ALL|wx.EXPAND,1)
3662    aSizer.Add(panel2,2, wx.ALL|wx.EXPAND,1)
3663    mainSizer.Add(aSizer,1, wx.ALL|wx.EXPAND,1)
3664
3665    # make OK button
3666    btnsizer = wx.BoxSizer(wx.HORIZONTAL)
3667    btn = wx.Button(dlg, wx.ID_CLOSE) 
3668    btn.Bind(wx.EVT_BUTTON,_onClose)
3669    btnsizer.Add(btn)
3670    mainSizer.Add(btnsizer, 0, wx.ALIGN_CENTER|wx.ALL, 5)
3671
3672    mainSizer.Fit(dlg)
3673    dlg.SetMinSize(dlg.GetSize())
3674    dlg.ShowModal()
3675    dlg.Destroy()
Note: See TracBrowser for help on using the repository browser.