source: trunk/GSASIIconstrGUI.py @ 4584

Last change on this file since 4584 was 4584, checked in by toby, 14 months ago

prevent duplicate names for rigid bodies; make sure vector RB has 3+ atoms (less causes error in Orient ref sel); fix res ComboBox? selection

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