source: trunk/GSASIIconstrGUI.py @ 4431

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

new code for locating RBs; redo GUI for RB extractor

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