source: trunk/GSASIIconstrGUI.py @ 3835

Last change on this file since 3835 was 3835, checked in by vondreele, 3 years ago

Allow constraints on modulation vector
changes to incommensurate mag. str factr math - a lot closer to correct

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 103.4 KB
Line 
1# -*- coding: utf-8 -*-
2#GSASIIconstrGUI - constraint GUI routines
3########### SVN repository information ###################
4# $Date: 2019-03-05 20:57:21 +0000 (Tue, 05 Mar 2019) $
5# $Author: vondreele $
6# $Revision: 3835 $
7# $URL: trunk/GSASIIconstrGUI.py $
8# $Id: GSASIIconstrGUI.py 3835 2019-03-05 20:57:21Z vondreele $
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 sys
19import wx
20import wx.grid as wg
21import random as ran
22import numpy as np
23import numpy.ma as ma
24import numpy.linalg as nl
25import os.path
26import GSASIIpath
27GSASIIpath.SetVersionNumber("$Revision: 3835 $")
28import GSASIIElem as G2elem
29import GSASIIElemGUI as G2elemGUI
30import GSASIIstrIO as G2stIO
31import GSASIImapvars as G2mv
32import GSASIImath as G2mth
33import GSASIIlattice as G2lat
34import GSASIIdataGUI as G2gd
35import GSASIIctrlGUI as G2G
36import GSASIIplot as G2plt
37import GSASIIobj as G2obj
38import GSASIIspc as G2spc
39VERY_LIGHT_GREY = wx.Colour(235,235,235)
40
41class ConstraintDialog(wx.Dialog):
42    '''Window to edit Constraint values
43    '''
44    def __init__(self,parent,title,text,data,separator='*',varname="",varyflag=False):
45        wx.Dialog.__init__(self,parent,-1,'Edit '+title, 
46            pos=wx.DefaultPosition,style=wx.DEFAULT_DIALOG_STYLE)
47        self.data = data[:]
48        self.newvar = [varname,varyflag]
49        panel = wx.Panel(self)
50        mainSizer = wx.BoxSizer(wx.VERTICAL)
51        topLabl = wx.StaticText(panel,-1,text)
52        mainSizer.Add((10,10),1)
53        mainSizer.Add(topLabl,0,wx.ALIGN_CENTER_VERTICAL|wx.LEFT,10)
54        mainSizer.Add((10,10),1)
55        dataGridSizer = wx.FlexGridSizer(cols=3,hgap=2,vgap=2)
56        self.OkBtn = wx.Button(panel,wx.ID_OK)
57        for id in range(len(self.data)):
58            lbl1 = lbl = str(self.data[id][1])
59            if lbl[-1] != '=': lbl1 = lbl + ' ' + separator + ' '
60            name = wx.StaticText(panel,wx.ID_ANY,lbl1,
61                                 style=wx.ALIGN_RIGHT)
62            scale = G2G.ValidatedTxtCtrl(panel,self.data[id],0,
63                                          typeHint=float,
64                                          OKcontrol=self.DisableOK)
65            dataGridSizer.Add(name,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,5)
66            dataGridSizer.Add(scale,0,wx.RIGHT,3)
67            if ':' in lbl:
68                dataGridSizer.Add(
69                    wx.StaticText(panel,-1,G2obj.fmtVarDescr(lbl)),
70                    0,wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,3)
71            else:
72                dataGridSizer.Add((-1,-1))
73        if title == 'New Variable':
74            name = wx.StaticText(panel,wx.ID_ANY,"New variable's\nname (optional)",
75                                 style=wx.ALIGN_CENTER)
76            scale = G2G.ValidatedTxtCtrl(panel,self.newvar,0,
77                                          typeHint=str,notBlank=False)
78            dataGridSizer.Add(name,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,5)
79            dataGridSizer.Add(scale,0,wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,3)
80            self.refine = wx.CheckBox(panel,label='Refine?')
81            self.refine.SetValue(self.newvar[1]==True)
82            self.refine.Bind(wx.EVT_CHECKBOX, self.OnCheckBox)
83            dataGridSizer.Add(self.refine,0,wx.RIGHT|wx.ALIGN_CENTER_VERTICAL,3)
84           
85        mainSizer.Add(dataGridSizer,0,wx.EXPAND)
86        self.OkBtn.Bind(wx.EVT_BUTTON, self.OnOk)
87        self.OkBtn.SetDefault()
88        cancelBtn = wx.Button(panel,wx.ID_CANCEL)
89        cancelBtn.Bind(wx.EVT_BUTTON, self.OnCancel)
90        btnSizer = wx.BoxSizer(wx.HORIZONTAL)
91        btnSizer.Add((20,20),1)
92        btnSizer.Add(self.OkBtn)
93        btnSizer.Add((20,20),1)
94        btnSizer.Add(cancelBtn)
95        btnSizer.Add((20,20),1)
96
97        mainSizer.Add(btnSizer,0,wx.EXPAND|wx.BOTTOM|wx.TOP, 10)
98        panel.SetSizer(mainSizer)
99        panel.Fit()
100        self.Fit()
101        self.CenterOnParent()
102       
103    def DisableOK(self,setting):
104        if setting:
105            self.OkBtn.Enable()
106        else:
107            self.OkBtn.Disable()
108
109    def OnCheckBox(self,event):
110        self.newvar[1] = self.refine.GetValue()
111
112    def OnOk(self,event):
113        parent = self.GetParent()
114        parent.Raise()
115        self.EndModal(wx.ID_OK)             
116
117    def OnCancel(self,event):
118        parent = self.GetParent()
119        parent.Raise()
120        self.EndModal(wx.ID_CANCEL)             
121
122    def GetData(self):
123        return self.data
124       
125################################################################################
126#####  Constraints
127################################################################################           
128       
129def UpdateConstraints(G2frame,data):
130    '''Called when Constraints tree item is selected.
131    Displays the constraints in the data window
132    '''
133       
134    def FindEquivVarb(name,nameList):
135        'Creates a list of variables appropriate to constrain with name'
136        outList = []
137        #phlist = []
138        items = name.split(':')
139        namelist = [items[2],]
140        if 'dA' in name:
141            namelist = ['dAx','dAy','dAz']
142        elif 'AU' in name:
143            namelist = ['AUiso','AU11','AU22','AU33','AU12','AU13','AU23']
144        elif 'AM' in name:
145            namelist = ['AMx','AMy','AMz']
146        elif items[-1] in ['A0','A1','A2','A3','A4','A5']:
147            namelist = ['A0','A1','A2','A3','A4','A5']
148        elif items[-1] in ['D11','D22','D33','D12','D13','D23']:
149            namelist = ['D11','D22','D33','D12','D13','D23']
150        elif 'Tm' in name:
151            namelist = ['Tmin','Tmax']
152        elif 'MX' in name or 'MY' in name or 'MZ' in name:
153            namelist = ['MXcos','MYcos','MZcos','MXsin','MYsin','MZsin']
154        elif 'mV' in name:
155            namelist = ['mV0','mV1','mV2']
156        elif 'RB' in name:
157            rbfx = 'RB'+items[2][2]
158            if 'T' in name and 'Tr' not in name:
159                namelist = [rbfx+'T11',rbfx+'T22',rbfx+'T33',rbfx+'T12',rbfx+'T13',rbfx+'T23']
160            if 'L' in name:
161                namelist = [rbfx+'L11',rbfx+'L22',rbfx+'L33',rbfx+'L12',rbfx+'L13',rbfx+'L23']
162            if 'S' in name:
163                namelist = [rbfx+'S12',rbfx+'S13',rbfx+'S21',rbfx+'S23',rbfx+'S31',rbfx+'S32',rbfx+'SAA',rbfx+'SBB']
164            if 'U' in name:
165                namelist = [rbfx+'U',]
166
167        for item in nameList:
168            keys = item.split(':')
169            #if keys[0] not in phlist:
170            #    phlist.append(keys[0])
171            if items[1] == '*' and keys[2] in namelist: # wildcard -- select only sequential options
172                keys[1] = '*'
173                mitem = ':'.join(keys)
174                if mitem == name: continue
175                if mitem not in outList: outList.append(mitem)
176            elif keys[2] in namelist and item != name:
177                outList.append(item)
178        return outList
179       
180    def SelectVarbs(page,FrstVarb,varList,legend,constType):
181        '''Select variables used in constraints after one variable has
182        been selected. This routine determines the appropriate variables to be
183        used based on the one that has been selected and asks for more to be added.
184
185        It then creates the constraint and adds it to the constraints list.
186       
187        Called from OnAddEquivalence, OnAddFunction & OnAddConstraint (all but
188        OnAddHold)
189
190        :param list page: defines calling page -- type of variables to be used
191        :parm GSASIIobj.G2VarObj FrstVarb: reference to first selected variable
192        :param list varList: list of other appropriate variables to select from
193        :param str legend: header for selection dialog
194        :param str constType: type of constraint to be generated
195        :returns: a constraint, as defined in
196          :ref:`GSASIIobj <Constraint_definitions_table>`
197        '''
198        choices = [[i]+list(G2obj.VarDescr(i)) for i in varList]
199        meaning = G2obj.getDescr(FrstVarb.name)
200        if not meaning:
201            meaning = "(no definition found!)"
202        l = str(FrstVarb).split(':')
203        # make lists of phases & histograms to iterate over
204        phaselist = [l[0]]
205        if l[0]:
206            phaselbl = ['phase #'+l[0]]
207            if len(Phases) > 1:
208                phaselist += ['all'] 
209                phaselbl += ['all phases']
210        else:
211            phaselbl = ['']
212        histlist = [l[1]]
213        if l[1] == '*':
214            pass
215        elif l[1]:
216            histlbl = ['histogram #'+l[1]]
217            if len(Histograms) > 1:
218                histlist += ['all']
219                histlbl += ['all histograms']
220                typ = Histograms[G2obj.LookupHistName(l[1])[0]]['Instrument Parameters'][0]['Type'][1]
221                i = 0
222                for hist in Histograms:
223                    if Histograms[hist]['Instrument Parameters'][0]['Type'][1] == typ: i += 1
224                if i > 1:
225                    histlist += ['all='+typ]
226                    histlbl += ['all '+typ+' histograms']
227        else:
228            histlbl = ['']
229        # make a list of equivalent parameter names
230        nameList = [FrstVarb.name]
231        for var in varList:
232            nam = var.split(":")[2]
233            if nam not in nameList: nameList += [nam]
234        # add "wild-card" names to the list of variables
235        if l[1] == '*':
236            pass
237        elif page[1] == 'phs':
238            if 'RB' in FrstVarb.name:
239                pass
240            elif FrstVarb.atom is None:
241                for nam in nameList:
242                    for ph,plbl in zip(phaselist,phaselbl):
243                        if plbl: plbl = 'For ' + plbl
244                        var = ph+"::"+nam
245                        if var == str(FrstVarb) or var in varList: continue
246                        varList += [var]
247                        choices.append([var,plbl,meaning])
248            else:
249                for nam in nameList:
250                    for ph,plbl in zip(phaselist,phaselbl):
251                        if plbl: plbl = ' in ' + plbl
252                        for atype in ['']+TypeList:
253                            if atype:
254                                albl = "For "+atype+" atoms"
255                                akey = "all="+atype                       
256                            else:
257                                albl = "For all atoms"
258                                akey = "all"
259                            var = ph+"::"+nam+":"+akey
260                            if var == str(FrstVarb) or var in varList: continue
261                            varList += [var]
262                            choices.append([var,albl+plbl,meaning])
263        elif page[1] == 'hap':
264            if FrstVarb.name == "Scale":
265                meaning = "Phase fraction"
266            for nam in nameList:
267                for ph,plbl in zip(phaselist,phaselbl):
268                    if plbl: plbl = 'For ' + plbl
269                    for hst,hlbl in zip(histlist,histlbl):
270                        if hlbl:
271                            if plbl:
272                                hlbl = ' in ' + hlbl
273                            else:
274                                hlbl = 'For ' + hlbl                               
275                            var = ph+":"+hst+":"+nam
276                            if var == str(FrstVarb) or var in varList: continue
277                            varList += [var]
278                            choices.append([var,plbl+hlbl,meaning])
279        elif page[1] == 'hst':
280            if FrstVarb.name == "Scale":
281                meaning = "Scale factor"
282            for nam in nameList:
283                for hst,hlbl in zip(histlist,histlbl):
284                    if hlbl:
285                        hlbl = 'For ' + hlbl                               
286                        var = ":"+hst+":"+nam
287                        if var == str(FrstVarb) or var in varList: continue
288                        varList += [var]
289                        choices.append([var,hlbl,meaning])
290        elif page[1] == 'glb' or page[1] == 'sym':
291            pass
292        else:
293            raise Exception('Unknown constraint page '+ page[1])                   
294        if len(choices):
295            l1 = l2 = 1
296            for i1,i2,i3 in choices:
297                l1 = max(l1,len(i1))
298                l2 = max(l2,len(i2))
299            fmt = "{:"+str(l1)+"s} {:"+str(l2)+"s} {:s}"
300            atchoices = [fmt.format(*i1) for i1 in choices] # reformat list as str with columns
301            dlg = G2G.G2MultiChoiceDialog(
302                G2frame,legend,
303                'Constrain '+str(FrstVarb)+' with...',atchoices,
304                toggle=False,size=(625,400),monoFont=True)
305            dlg.CenterOnParent()
306            res = dlg.ShowModal()
307            Selections = dlg.GetSelections()[:]
308            dlg.Destroy()
309            if res != wx.ID_OK: return []
310            if len(Selections) == 0:
311                dlg = wx.MessageDialog(
312                    G2frame,
313                    'No variables were selected to include with '+str(FrstVarb),
314                    'No variables')
315                dlg.CenterOnParent()
316                dlg.ShowModal()
317                dlg.Destroy()
318                return []
319        else:
320            dlg = wx.MessageDialog(
321                G2frame,
322                'There are no appropriate variables to include with '+str(FrstVarb),
323                'No variables')
324            dlg.CenterOnParent()
325            dlg.ShowModal()
326            dlg.Destroy()
327            return []
328        # now process the variables provided by the user
329        varbs = [str(FrstVarb),] # list of selected variables
330        for sel in Selections:
331            var = varList[sel]
332            # phase(s) included
333            l = var.split(':')
334            if l[0] == "all":
335                phlist = [str(Phases[phase]['pId']) for phase in Phases]
336            else:
337                phlist = [l[0]]
338            # histogram(s) included
339            if l[1] == "all":
340                hstlist = [str(Histograms[hist]['hId']) for hist in Histograms]
341            elif '=' in l[1]:
342                htyp = l[1].split('=')[1]
343                hstlist = [str(Histograms[hist]['hId']) for hist in Histograms if 
344                           Histograms[hist]['Instrument Parameters'][0]['Type'][1] == htyp]
345            else:
346                hstlist = [l[1]]
347            if len(l) == 3:
348                for ph in phlist:
349                    for hst in hstlist:
350                        var = ph + ":" + hst + ":" + l[2]
351                        if var in varbs: continue
352                        varbs.append(var)
353            else: # constraints with atoms or rigid bodies
354                if len(l) == 5: # rigid body parameter
355                    var = ':'.join(l)
356                    if var in varbs: continue
357                    varbs.append(var)
358                elif l[3] == "all":
359                    for ph in phlist:
360                        key = G2obj.LookupPhaseName(ph)[0]
361                        for hst in hstlist: # should be blank
362                            for iatm,at in enumerate(Phases[key]['Atoms']):
363                                var = ph + ":" + hst + ":" + l[2] + ":" + str(iatm)
364                                if var in varbs: continue
365                                varbs.append(var)
366                elif '=' in l[3]:
367                    for ph in phlist:
368                        key = G2obj.LookupPhaseName(ph)[0]
369                        cx,ct,cs,cia = Phases[key]['General']['AtomPtrs']
370                        for hst in hstlist: # should be blank
371                            atyp = l[3].split('=')[1]
372                            for iatm,at in enumerate(Phases[key]['Atoms']):
373                                if at[ct] != atyp: continue
374                                var = ph + ":" + hst + ":" + l[2] + ":" + str(iatm)
375                                if var in varbs: continue
376                                varbs.append(var)
377                else:
378                    for ph in phlist:
379                        key = G2obj.LookupPhaseName(ph)[0]
380                        for hst in hstlist: # should be blank
381                            var = ph + ":" + hst + ":" + l[2] + ":" + l[3]
382                            if var in varbs: continue
383                            varbs.append(var)
384        if len(varbs) >= 1 or 'constraint' in constType:
385            constr = [[1.0,FrstVarb]]
386            for item in varbs[1:]:
387                constr += [[1.0,G2obj.G2VarObj(item)]]
388            if 'equivalence' in constType:
389                return [constr+[None,None,'e']]
390            elif 'function' in constType:
391                return [constr+[None,False,'f']]
392            elif 'constraint' in constType:
393                return [constr+[1.0,None,'c']]
394            else:
395                raise Exception('Unknown constraint type: '+str(constType))
396        else:
397            dlg = wx.MessageDialog(
398                G2frame,
399                'There are no selected variables to include with '+str(FrstVarb),
400                'No variables')
401            dlg.CenterOnParent()
402            dlg.ShowModal()
403            dlg.Destroy()
404        return []
405   
406    def ConstraintsLoad(data,newcons=[]):
407        '''Load all constraints. Constraints based on symmetry (etc.)
408        are generated by running :func:`GSASIIstrIO.GetPhaseData`.
409        '''
410        G2mv.InitVars()
411        #Find all constraints
412        constraintSet = []
413        for key in data:
414            if key.startswith('_'): continue
415            constraintSet += data[key]
416        if newcons:
417            constraintSet = constraintSet + newcons
418        constDictList,fixedList,ignored = G2stIO.ProcessConstraints(constraintSet)
419        # generate symmetry constraints to check for conflicts
420        rigidbodyDict = G2frame.GPXtree.GetItemPyData(   
421            G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Rigid bodies'))
422        rbIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[]})
423        rbVary,rbDict = G2stIO.GetRigidBodyModels(rigidbodyDict,Print=False)
424        Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables,MFtables,maxSSwave = G2stIO.GetPhaseData(
425            Phases,RestraintDict=None,rbIds=rbIds,Print=False) # generates atom symmetry constraints
426        return constDictList,phaseDict,fixedList
427           
428    def ConstraintsCheck(data,newcons=[]):
429        '''Load constraints & check them for errors. Since error checking
430        can cause changes in constraints in case of repairable conflicts
431        between equivalences, reload the constraints again after the check.
432        This could probably be done more effectively, but only reloading when
433        needed, but the reload should be quick.
434        '''
435        constDictList,phaseDict,fixedList = ConstraintsLoad(data,newcons)
436        msg = G2mv.EvaluateMultipliers(constDictList,phaseDict)
437        if msg:
438            return 'Unable to interpret multiplier(s): '+msg,''
439        res = G2mv.CheckConstraints('',constDictList,fixedList)
440        # reload constraints in case any were merged in MoveConfEquiv
441        ConstraintsLoad(data,newcons)
442        return res
443
444    def CheckAddedConstraint(newcons):
445        '''Check a new constraint that has just been input.
446        If there is an error display a message and discard the last entry
447
448        Since the varylist is not available, no warning messages
449        should be generated here
450
451        :returns: True if constraint should be added
452        '''
453       
454        errmsg,warnmsg = ConstraintsCheck(data,newcons)
455        if errmsg:
456            G2frame.ErrorDialog('Constraint Error',
457                'Error with newly added constraint:\n'+errmsg+
458                '\nIgnoring newly added constraint',parent=G2frame)
459            # reset error status
460            errmsg,warnmsg = ConstraintsCheck(data)
461            if errmsg:
462                print (errmsg)
463                print (G2mv.VarRemapShow([],True))
464            return False
465        elif warnmsg:
466            print ('Unexpected contraint warning:\n'+warnmsg)
467        return True
468
469    def CheckChangedConstraint():
470        '''Check all constraints after an edit has been made.
471        If there is an error display a message and reject the change.
472
473        Since the varylist is not available, no warning messages
474        should be generated.
475       
476        :returns: True if the edit should be retained
477        '''
478        errmsg,warnmsg = ConstraintsCheck(data)
479        if errmsg:
480            G2frame.ErrorDialog('Constraint Error',
481                'Error after editing constraint:\n'+errmsg+
482                '\nDiscarding last constraint edit',parent=G2frame)
483            # reset error status
484            errmsg,warnmsg = ConstraintsCheck(data)
485            if errmsg:
486                print (errmsg)
487                print (G2mv.VarRemapShow([],True))
488            return False
489        elif warnmsg:
490            print ('Unexpected contraint warning:\n'+warnmsg)
491        return True
492             
493    def PageSelection(page):
494        'Decode page reference'
495        if page[1] == "phs":
496            vartype = "phase"
497            varList = G2obj.removeNonRefined(phaseList)  # remove any non-refinable prms from list
498            constrDictEnt = 'Phase'
499        elif page[1] == "hap":
500            vartype = "Histogram*Phase"
501            varList = G2obj.removeNonRefined(hapList)  # remove any non-refinable prms from list
502            constrDictEnt = 'HAP'
503        elif page[1] == "hst":
504            vartype = "Histogram"
505            varList = G2obj.removeNonRefined(histList)  # remove any non-refinable prms from list
506            constrDictEnt = 'Hist'
507        elif page[1] == "glb":
508            vartype = "Global"
509            varList = G2obj.removeNonRefined(globalList)   # remove any non-refinable prms from list
510
511            constrDictEnt = 'Global'
512        elif page[1] == "sym":
513            return None,None,None
514        else:
515            raise Exception('Should not happen!')
516        return vartype,varList,constrDictEnt
517
518    def OnAddHold(event):
519        '''Create a new Hold constraint
520
521        Hold constraints allows the user to select one variable (the list of available
522        variables depends on which tab is currently active).
523        '''
524        page = G2frame.Page
525        vartype,varList,constrDictEnt = PageSelection(page)
526        if vartype is None: return
527        varList = G2obj.SortVariables(varList)
528        title1 = "Hold "+vartype+" variable"
529        if not varList:
530            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
531                parent=G2frame)
532            return
533        l2 = l1 = 1
534        for i in varList:
535            l1 = max(l1,len(i))
536            loc,desc = G2obj.VarDescr(i)
537            l2 = max(l2,len(loc))
538        fmt = "{:"+str(l1)+"s} {:"+str(l2)+"s} {:s}"
539        varListlbl = [fmt.format(i,*G2obj.VarDescr(i)) for i in varList]
540        #varListlbl = ["("+i+") "+G2obj.fmtVarDescr(i) for i in varList]
541        legend = "Select variables to hold (Will not be varied, even if vary flag is set)"
542        dlg = G2G.G2MultiChoiceDialog(G2frame,
543            legend,title1,varListlbl,toggle=False,size=(625,400),monoFont=True)
544        dlg.CenterOnParent()
545        if dlg.ShowModal() == wx.ID_OK:
546            for sel in dlg.GetSelections():
547                Varb = varList[sel]
548                VarObj = G2obj.G2VarObj(Varb)
549                newcons = [[[0.0,VarObj],None,None,'h']]
550                if CheckAddedConstraint(newcons):
551                    data[constrDictEnt] += newcons
552        dlg.Destroy()
553        wx.CallAfter(OnPageChanged,None)
554       
555    def OnAddEquivalence(event):
556        '''add an Equivalence constraint'''
557        page = G2frame.Page
558        vartype,varList,constrDictEnt = PageSelection(page)
559        if vartype is None: return
560        title1 = "Create equivalence constraint between "+vartype+" variables"
561        title2 = "Select additional "+vartype+" variable(s) to be equivalent with "
562        if not varList:
563            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
564                parent=G2frame)
565            return
566#        legend = "Select variables to make equivalent (only one of the variables will be varied when all are set to be varied)"
567        GetAddVars(page,title1,title2,varList,constrDictEnt,'equivalence')
568       
569    def OnAddAtomEquiv(event):
570        ''' Add equivalences between all parameters on atoms '''
571        page = G2frame.Page
572        vartype,varList,constrDictEnt = PageSelection(page)
573        if vartype is None: return
574        title1 = "Setup equivalent atom variables"
575        title2 = "Select additional atoms(s) to be equivalent with "
576        if not varList:
577            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
578                parent=G2frame)
579            return
580#        legend = "Select atoms to make equivalent (only one of the atom variables will be varied when all are set to be varied)"
581        GetAddAtomVars(page,title1,title2,varList,constrDictEnt,'equivalence')
582       
583    def OnAddRiding(event):
584        ''' Add riding equivalences between all parameters on atoms '''
585        page = G2frame.Page
586        vartype,varList,constrDictEnt = PageSelection(page)
587        if vartype is None: return
588        title1 = "Setup riding atoms "
589        title2 = "Select additional atoms(s) to ride on "
590        if not varList:
591            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
592                parent=G2frame)
593            return
594#        legend = "Select atoms to ride (only one of the atom variables will be varied when all are set to be varied)"
595        GetAddAtomVars(page,title1,title2,varList,constrDictEnt,'riding')
596   
597    def OnAddFunction(event):
598        '''add a Function (new variable) constraint'''
599        page = G2frame.Page
600        vartype,varList,constrDictEnt = PageSelection(page)
601        if vartype is None: return
602        title1 = "Setup new variable based on "+vartype+" variables"
603        title2 = "Include additional "+vartype+" variable(s) to be included with "
604        if not varList:
605            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
606                parent=G2frame)
607            return
608#        legend = "Select variables to include in a new variable (the new variable will be varied when all included variables are varied)"
609        GetAddVars(page,title1,title2,varList,constrDictEnt,'function')
610                       
611    def OnAddConstraint(event):
612        '''add a constraint equation to the constraints list'''
613        page = G2frame.Page
614        vartype,varList,constrDictEnt = PageSelection(page)
615        if vartype is None: return
616        title1 = "Creating constraint on "+vartype+" variables"
617        title2 = "Select additional "+vartype+" variable(s) to include in constraint with "
618        if not varList:
619            G2frame.ErrorDialog('No variables','There are no variables of type '+vartype,
620                parent=G2frame)
621            return
622#        legend = "Select variables to include in a constraint equation (the values will be constrainted to equal a specified constant)"
623        GetAddVars(page,title1,title2,varList,constrDictEnt,'constraint')
624
625    def GetAddVars(page,title1,title2,varList,constrDictEnt,constType):
626        '''Get the variables to be added for OnAddEquivalence, OnAddFunction,
627        and OnAddConstraint. Then create and check the constraint.
628        '''
629        #varListlbl = ["("+i+") "+G2obj.fmtVarDescr(i) for i in varList]
630        if constType == 'equivalence':
631            omitVars = G2mv.GetDependentVars()
632        else:
633            omitVars = []
634        varList = G2obj.SortVariables([i for i in varList if i not in omitVars])
635        l2 = l1 = 1
636        for i in varList:
637            l1 = max(l1,len(i))
638            loc,desc = G2obj.VarDescr(i)
639            l2 = max(l2,len(loc))
640        fmt = "{:"+str(l1)+"s} {:"+str(l2)+"s} {:s}"
641        varListlbl = [fmt.format(i,*G2obj.VarDescr(i)) for i in varList]
642        dlg = G2G.G2SingleChoiceDialog(G2frame,'Select 1st variable:',
643            title1,varListlbl,monoFont=True,size=(625,400))
644        dlg.CenterOnParent()
645        if dlg.ShowModal() == wx.ID_OK:
646            if constType == 'equivalence':
647                omitVars = G2mv.GetDependentVars() + G2mv.GetIndependentVars()
648            sel = dlg.GetSelection()
649            FrstVarb = varList[sel]
650            VarObj = G2obj.G2VarObj(FrstVarb)
651            moreVarb = G2obj.SortVariables(FindEquivVarb(FrstVarb,[i for i in varList if i not in omitVars]))
652            newcons = SelectVarbs(page,VarObj,moreVarb,title2+FrstVarb,constType)
653            if len(newcons) > 0:
654                if CheckAddedConstraint(newcons):
655                    data[constrDictEnt] += newcons
656        dlg.Destroy()
657        wx.CallAfter(OnPageChanged,None)
658                       
659    def FindNeighbors(phase,FrstName,AtNames):
660        General = phase['General']
661        cx,ct,cs,cia = General['AtomPtrs']
662        Atoms = phase['Atoms']
663        atNames = [atom[ct-1] for atom in Atoms]
664        Cell = General['Cell'][1:7]
665        Amat,Bmat = G2lat.cell2AB(Cell)
666        atTypes = General['AtomTypes']
667        Radii = np.array(General['BondRadii'])
668        AtInfo = dict(zip(atTypes,Radii)) #or General['BondRadii']
669        Orig = atNames.index(FrstName.split()[1])
670        OType = Atoms[Orig][ct]
671        XYZ = G2mth.getAtomXYZ(Atoms,cx)       
672        Neigh = []
673        Dx = np.inner(Amat,XYZ-XYZ[Orig]).T
674        dist = np.sqrt(np.sum(Dx**2,axis=1))
675        sumR = AtInfo[OType]+0.5    #H-atoms only!
676        IndB = ma.nonzero(ma.masked_greater(dist-0.85*sumR,0.))
677        for j in IndB[0]:
678            if j != Orig:
679                Neigh.append(AtNames[j])
680        return Neigh
681       
682    def GetAddAtomVars(page,title1,title2,varList,constrDictEnt,constType):
683        '''Get the atom variables to be added for OnAddAtomEquiv. Then create and
684        check the constraints. Riding for H atoms only.
685        '''
686        Atoms = {G2obj.VarDescr(i)[0]:[] for i in varList if 'Atom' in G2obj.VarDescr(i)[0]}
687        for item in varList:
688            atName = G2obj.VarDescr(item)[0]
689            if atName in Atoms:
690                Atoms[atName].append(item)
691        AtNames = list(Atoms.keys())
692        AtNames.sort()
693        dlg = G2G.G2SingleChoiceDialog(G2frame,'Select 1st atom:',
694            title1,AtNames,monoFont=True,size=(625,400))
695        dlg.CenterOnParent()
696        FrstAtom = ''
697        if dlg.ShowModal() == wx.ID_OK:
698            sel = dlg.GetSelection()
699            FrstAtom = AtNames[sel]
700            if 'riding' in constType:
701                phaseName = (FrstAtom.split(' in ')[1]).strip()
702                phase = Phases[phaseName]
703                AtNames = FindNeighbors(phase,FrstAtom,AtNames)
704            else:
705                AtNames.remove(FrstAtom)
706        dlg.Destroy()
707        if FrstAtom == '':
708            print ('no atom selected')
709            return
710        dlg = G2G.G2MultiChoiceDialog(
711            G2frame,title2+FrstAtom,
712            'Constrain '+str(FrstAtom)+' with...',AtNames,
713            toggle=False,size=(625,400),monoFont=True)
714        if dlg.ShowModal() == wx.ID_OK:
715            Selections = dlg.GetSelections()[:]
716        else:
717            print ('no target atom selected')
718            dlg.Destroy()
719            return
720        dlg.Destroy()
721        for name in Atoms[FrstAtom]:
722            newcons = []
723            constr = []
724            if 'riding' in constType:
725                if 'AUiso' in name:
726                    constr = [[1.0,G2obj.G2VarObj(name)]]
727                elif 'AU11' in name:
728                    pass
729                elif 'AU' not in name:
730                    constr = [[1.0,G2obj.G2VarObj(name)]]
731            else:
732                constr = [[1.0,G2obj.G2VarObj(name)]]
733            pref = ':'+name.rsplit(':',1)[0].split(':',1)[1]    #get stuff between phase id & atom id
734            for sel in Selections:
735                name2 = Atoms[AtNames[sel]][0]
736                pid = name2.split(':',1)[0]                     #get phase id for 2nd atom
737                id = name2.rsplit(':',1)[-1]                    #get atom no. for 2nd atom
738                if 'riding' in constType:
739                    pref = pid+pref
740                    if 'AUiso' in pref:
741                        parts = pref.split('AUiso')
742                        constr += [[1.2,G2obj.G2VarObj('%s:%s'%(parts[0]+'AUiso',id))]]
743                    elif 'AU' not in pref:
744                        constr += [[1.0,G2obj.G2VarObj('%s:%s'%(pref,id))]]
745                else:
746                    constr += [[1.0,G2obj.G2VarObj('%s:%s'%(pid+pref,id))]]
747            if not constr:
748                continue
749            if 'frac' in pref and 'riding' not in constType:
750                newcons = [constr+[1.0,None,'c']]
751            else:
752                newcons = [constr+[None,None,'e']]
753            if len(newcons) > 0:
754                if CheckAddedConstraint(newcons):
755                    data[constrDictEnt] += newcons
756        wx.CallAfter(OnPageChanged,None)
757                       
758    def MakeConstraintsSizer(name,pageDisplay):
759        '''Creates a sizer displaying all of the constraints entered of
760        the specified type.
761
762        :param str name: the type of constraints to be displayed ('HAP',
763          'Hist', 'Phase', 'Global', 'Sym-Generated')
764        :param wx.Panel pageDisplay: parent panel for sizer
765        :returns: wx.Sizer created by method
766        '''
767        if name == 'Sym-Generated':         #show symmetry generated constraints
768            Sizer1 =  wx.BoxSizer(wx.VERTICAL)
769            Sizer1.Add(wx.StaticText(pageDisplay,wx.ID_ANY,
770                                    'Equivalences generated based on cell/space group input'))
771            Sizer1.Add((-1,5))
772            Sizer = wx.FlexGridSizer(0,2,0,0)
773            Sizer1.Add(Sizer)
774            for sym in G2mv.GetSymEquiv():
775                Sizer.Add(wx.StaticText(pageDisplay,wx.ID_ANY,'EQUIV'),
776                           0,wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER|wx.RIGHT|wx.LEFT,3)
777                Sizer.Add(wx.StaticText(pageDisplay,wx.ID_ANY,sym))
778                Sizer.Add((-1,-1))
779                Sizer.Add((-1,2))
780            return Sizer1
781        constSizer = wx.FlexGridSizer(0,6,0,0)
782        maxlen = 50 # characters before wrapping a constraint
783        for Id,item in enumerate(data[name]):
784            refineflag = False
785            helptext = ""
786            eqString = ['',]
787            problemItem = False
788            for term in item[:-3]:
789                if str(term[1]) in G2mv.problemVars:
790                    problemItem = True
791            if item[-1] == 'h': # Hold on variable
792                constSizer.Add((-1,-1),0)              # blank space for edit button
793                typeString = 'FIXED'
794                var = str(item[0][1])
795                varMean = G2obj.fmtVarDescr(var)
796                eqString[-1] =  var +'   '
797                helptext = "Prevents variable:\n"+ var + " ("+ varMean + ")\nfrom being changed"
798            elif item[-1] == 'f' or item[-1] == 'e' or item[-1] == 'c': # not true on original-style (2011?) constraints
799                constEdit = wx.Button(pageDisplay,wx.ID_ANY,'Edit',style=wx.BU_EXACTFIT)
800                constEdit.Bind(wx.EVT_BUTTON,OnConstEdit)
801                constSizer.Add(constEdit,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)            # edit button
802                Indx[constEdit.GetId()] = [Id,name]
803                if item[-1] == 'f':
804                    helptext = "A new variable"
805                    if item[-3]:
806                        helptext += " named "+str(item[-3])
807                    helptext += " is created from a linear combination of the following variables:\n"
808                    for term in item[:-3]:
809                        var = str(term[1])
810                        if len(eqString[-1]) > maxlen:
811                            eqString.append(' ')
812                        m = term[0]
813                        if eqString[-1] != '':
814                            if m >= 0:
815                                eqString[-1] += ' + '
816                            else:
817                                eqString[-1] += ' - '
818                                m = abs(m)
819                        if m == 1:
820                            eqString[-1] += '{:} '.format(var)
821                        else:
822                            eqString[-1] += '{:g}*{:} '.format(m,var)
823                        varMean = G2obj.fmtVarDescr(var)
824                        helptext += "\n" + var + " ("+ varMean + ")"
825                    if '_Explain' in data:
826                        if data['_Explain'].get(item[-3]):
827                            helptext += '\n\n'
828                            helptext += data['_Explain'][item[-3]]
829                    # typeString = 'NEWVAR'
830                    # if item[-3]:
831                    #     eqString[-1] += ' = '+item[-3]
832                    # else:
833                    #     eqString[-1] += ' = New Variable'
834                    if item[-3]:
835                        typeString = item[-3] + ' = '
836                    else:
837                        typeString = 'New Variable = '
838                    #print 'refine',item[-2]
839                    refineflag = True
840                elif item[-1] == 'c':
841                    helptext = "The following variables constrained to equal a constant:"
842                    for term in item[:-3]:
843                        var = str(term[1])
844                        if len(eqString[-1]) > maxlen:
845                            eqString.append(' ')
846                        m = term[0]
847                        if eqString[-1] != '':
848                            if term[0] > 0:
849                                eqString[-1] += ' + '
850                            else:
851                                eqString[-1] += ' - '
852                                m = -term[0]
853                        if m == 1:
854                            eqString[-1] += '{:} '.format(var)
855                        else:
856                            eqString[-1] += '{:g}*{:} '.format(m,var)
857                        varMean = G2obj.fmtVarDescr(var)
858                        helptext += "\n" + var + " ("+ varMean + ")"
859                    typeString = 'CONST'
860                    eqString[-1] += ' = '+str(item[-3])
861                elif item[-1] == 'e':
862                    helptext = "The following variables are set to be equivalent, noting multipliers:"
863                    normval = item[:-3][1][0]
864                    for i,term in enumerate(item[:-3]):
865                        var = str(term[1])
866                        if term[0] == 0: term[0] = 1.0
867                        if len(eqString[-1]) > maxlen:
868                            eqString.append(' ')
869                        if i == 0: # move independent variable to end
870                            indepterm = term
871                            continue
872                        elif eqString[-1] != '':
873                            eqString[-1] += ' = '
874                        if normval/term[0] == 1:
875                            eqString[-1] += '{:}'.format(var)
876                        else:
877                            eqString[-1] += '{:g}*{:} '.format(normval/term[0],var)
878                        varMean = G2obj.fmtVarDescr(var)
879                        helptext += "\n" + var + " ("+ varMean + ")"
880                    if normval/indepterm[0] == 1:
881                        eqString[-1] += ' = {:} '.format(indepterm[1])
882                    else:
883                        eqString[-1] += ' = {:g}*{:} '.format(normval/indepterm[0],str(indepterm[1]))
884                    typeString = 'EQUIV'
885                else:
886                    print ('Unexpected constraint'+item)
887               
888            else:
889                print ('Removing old-style constraints')
890                data[name] = []
891                return constSizer
892            constDel = wx.Button(pageDisplay,wx.ID_ANY,'Delete',style=wx.BU_EXACTFIT)
893            constDel.Bind(wx.EVT_BUTTON,OnConstDel)
894            Indx[constDel.GetId()] = [Id,name]
895            constSizer.Add(constDel,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)             # delete button
896            if helptext:
897                ch = G2G.HelpButton(pageDisplay,helptext)
898                constSizer.Add(ch,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)
899            else:
900                constSizer.Add((-1,-1))
901            if refineflag:
902                ch = G2G.G2CheckBox(pageDisplay,'',item,-2)
903                constSizer.Add(ch,0,wx.LEFT|wx.RIGHT|wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER,1)
904            else:
905                constSizer.Add((-1,-1))               
906            constSizer.Add(wx.StaticText(pageDisplay,wx.ID_ANY,typeString),
907                           0,wx.ALIGN_CENTER_VERTICAL|wx.ALIGN_CENTER|wx.RIGHT|wx.LEFT,3)
908            if problemItem: eqString[-1] += ' -- Conflict: see console'
909            if len(eqString) > 1:
910                Eq = wx.BoxSizer(wx.VERTICAL)
911                for s in eqString:
912                    line = wx.StaticText(pageDisplay,wx.ID_ANY,s)
913                    if problemItem: line.SetBackgroundColour(wx.YELLOW)
914                    Eq.Add(line,0,wx.ALIGN_CENTER_VERTICAL)
915                Eq.Add((-1,4))
916            else:
917                Eq = wx.StaticText(pageDisplay,wx.ID_ANY,eqString[0])
918                if problemItem: Eq.SetBackgroundColour(wx.YELLOW)
919            constSizer.Add(Eq,1,wx.ALIGN_CENTER_VERTICAL)
920        return constSizer
921               
922    def OnConstDel(event):
923        'Delete a constraint'
924        Obj = event.GetEventObject()
925        Id,name = Indx[Obj.GetId()]
926        del data[name][Id]
927        ConstraintsLoad(data)
928        wx.CallAfter(OnPageChanged,None)
929
930    def OnConstEdit(event):
931        '''Called to edit an individual contraint in response to a
932        click on its Edit button
933        '''
934        Obj = event.GetEventObject()
935        Id,name = Indx[Obj.GetId()]
936        if data[name][Id][-1] == 'f':
937            items = data[name][Id][:-3]
938            constType = 'New Variable'
939            if data[name][Id][-3]:
940                varname = data[name][Id][-3]
941            else:
942                varname = ""
943            lbl = 'Enter value for each term in constraint; sum = new variable'
944            dlg = ConstraintDialog(G2frame,constType,lbl,items,
945                                   varname=varname,varyflag=data[name][Id][-2])
946        elif data[name][Id][-1] == 'c':
947            items = data[name][Id][:-3]+[
948                [data[name][Id][-3],'fixed value =']]
949            constType = 'Constraint'
950            lbl = 'Edit value for each term in constant constraint sum'
951            dlg = ConstraintDialog(G2frame,constType,lbl,items)
952        elif data[name][Id][-1] == 'e':
953            items = data[name][Id][:-3]
954            constType = 'Equivalence'
955            lbl = 'The following terms are set to be equal:'
956            dlg = ConstraintDialog(G2frame,constType,lbl,items,'/')
957        else:
958            return
959        try:
960            prev = data[name][Id][:]
961            if dlg.ShowModal() == wx.ID_OK:
962                result = dlg.GetData()
963                for i in range(len(data[name][Id][:-3])):
964                    if type(data[name][Id][i]) is tuple: # fix non-mutable construct
965                        data[name][Id][i] = list(data[name][Id][i])
966                    data[name][Id][i][0] = result[i][0]
967                if data[name][Id][-1] == 'c':
968                    data[name][Id][-3] = str(result[-1][0])
969                elif data[name][Id][-1] == 'f':
970                    # process the variable name to put in global form (::var)
971                    varname = str(dlg.newvar[0]).strip().replace(' ','_')
972                    if varname.startswith('::'):
973                        varname = varname[2:]
974                    varname = varname.replace(':',';')
975                    if varname:
976                        data[name][Id][-3] = varname
977                    else:
978                        data[name][Id][-3] = ''
979                    data[name][Id][-2] = dlg.newvar[1]
980                if not CheckChangedConstraint():
981                    data[name][Id] = prev
982        except:
983            import traceback
984            print (traceback.format_exc())
985        finally:
986            dlg.Destroy()
987        wx.CallAfter(OnPageChanged,None)
988   
989    def UpdateConstraintPanel(panel,typ):
990        '''Update the contents of the selected Constraint
991        notebook tab. Called in :func:`OnPageChanged`
992        '''
993        if panel.GetSizer(): panel.GetSizer().Clear(True)
994        Siz = wx.BoxSizer(wx.VERTICAL)
995        Siz.Add((5,5),0)
996        Siz.Add(MakeConstraintsSizer(typ,panel),1,wx.EXPAND)
997        panel.SetSizer(Siz,True)
998        Size = Siz.GetMinSize()
999        Size[0] += 40
1000        Size[1] = max(Size[1],450) + 20
1001        panel.SetSize(Size)
1002        panel.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
1003        panel.Show()
1004
1005    def OnPageChanged(event):
1006        '''Called when a tab is pressed or when a "select tab" menu button is
1007        used (see RaisePage), or to refresh the current tab contents (event=None)
1008        '''
1009        if event:       #page change event!
1010            page = event.GetSelection()
1011        else: # called directly, get current page
1012            page = G2frame.constr.GetSelection()
1013        G2frame.constr.ChangeSelection(page)
1014        text = G2frame.constr.GetPageText(page)
1015        G2frame.dataWindow.ConstraintEdit.Enable(G2G.wxID_EQUIVALANCEATOMS,False)
1016#        G2frame.dataWindow.ConstraintEdit.Enable(G2G.wxID_ADDRIDING,False)
1017        if text == 'Histogram/Phase':
1018            enableEditCons = [False]+4*[True]
1019            G2frame.Page = [page,'hap']
1020            UpdateConstraintPanel(HAPConstr,'HAP')
1021        elif text == 'Histogram':
1022            enableEditCons = [False]+4*[True]
1023            G2frame.Page = [page,'hst']
1024            UpdateConstraintPanel(HistConstr,'Hist')
1025        elif text == 'Phase':
1026            enableEditCons = 5*[True]
1027            G2frame.Page = [page,'phs']
1028            G2frame.dataWindow.ConstraintEdit.Enable(G2G.wxID_EQUIVALANCEATOMS,True)
1029#            G2frame.dataWindow.ConstraintEdit.Enable(G2G.wxID_ADDRIDING,True)
1030            if 'DELETED' in str(PhaseConstr):   #seems to be no other way to do this (wx bug)
1031                if GSASIIpath.GetConfigValue('debug'):
1032                    print ('DBG_wx error: PhaseConstr not cleanly deleted after Refine')
1033                return
1034            UpdateConstraintPanel(PhaseConstr,'Phase')
1035        elif text == 'Global':
1036            enableEditCons = [False]+4*[True]
1037            G2frame.Page = [page,'glb']
1038            UpdateConstraintPanel(GlobalConstr,'Global')
1039        else:
1040            enableEditCons = 5*[False]
1041            G2frame.Page = [page,'sym']
1042            UpdateConstraintPanel(SymConstr,'Sym-Generated')
1043        # remove menu items when not allowed
1044        for obj,flag in zip(G2frame.dataWindow.ConstraintEdit.GetMenuItems(),enableEditCons): 
1045            obj.Enable(flag)
1046        G2frame.dataWindow.SetDataSize()
1047
1048    def RaisePage(event):
1049        'Respond to a "select tab" menu button'
1050        try:
1051            i = (G2G.wxID_CONSPHASE,
1052                 G2G.wxID_CONSHAP,
1053                 G2G.wxID_CONSHIST,
1054                 G2G.wxID_CONSGLOBAL,
1055                 G2G.wxID_CONSSYM,
1056                ).index(event.GetId())
1057            G2frame.constr.SetSelection(i)
1058            wx.CallAfter(OnPageChanged,None)
1059        except ValueError:
1060            print('Unexpected event in RaisePage')
1061
1062    def SetStatusLine(text):
1063        G2frame.GetStatusBar().SetStatusText(text,1)
1064
1065    # UpdateConstraints execution starts here
1066    if not data:
1067        data.update({'Hist':[],'HAP':[],'Phase':[],'Global':[]})       #empty dict - fill it
1068    if 'Global' not in data:                                            #patch
1069        data['Global'] = []
1070    # DEBUG code ########################################
1071    #import GSASIIconstrGUI
1072    #reload(GSASIIconstrGUI)
1073    #reload(G2obj)
1074    #reload(G2stIO)
1075    #import GSASIIstrMain
1076    #reload(GSASIIstrMain)   
1077    #reload(G2mv)
1078    #reload(G2gd)
1079    ###################################################
1080    Histograms,Phases = G2frame.GetUsedHistogramsAndPhasesfromTree()
1081    if not len(Phases) or not len(Histograms):
1082        dlg = wx.MessageDialog(G2frame,'You need both phases and histograms to see Constraints',
1083            'No phases or histograms')
1084        dlg.CenterOnParent()
1085        dlg.ShowModal()
1086        dlg.Destroy()
1087        return
1088    G2obj.IndexAllIds(Histograms,Phases)
1089    ##################################################################################
1090    # patch: convert old-style (str) variables in constraints to G2VarObj objects
1091    for key,value in data.items():
1092        if key.startswith('_'): continue
1093        j = 0
1094        for cons in value:
1095            #print cons             # DEBUG
1096            for i in range(len(cons[:-3])):
1097                if type(cons[i][1]) is str:
1098                    cons[i][1] = G2obj.G2VarObj(cons[i][1])
1099                    j += 1
1100        if j:
1101            print (str(key) + ': '+str(j)+' variable(s) as strings converted to objects')
1102    ##################################################################################
1103    rigidbodyDict = G2frame.GPXtree.GetItemPyData(
1104        G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Rigid bodies'))
1105    rbIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[]})
1106    rbVary,rbDict = G2stIO.GetRigidBodyModels(rigidbodyDict,Print=False)
1107    badPhaseParms = ['Ax','Ay','Az','Amul','AI/A','Atype','SHorder','AwaveType','FwaveType','PwaveType','MwaveType','Vol','isMag',]
1108    globalList = list(rbDict.keys())
1109    globalList.sort()
1110    try:
1111        AtomDict = dict([Phases[phase]['pId'],Phases[phase]['Atoms']] for phase in Phases)
1112    except KeyError:
1113        G2frame.ErrorDialog('Constraint Error','Constraints cannot be set until a cycle of least squares'+
1114                            ' has been run.\nWe suggest you refine a scale factor.')
1115        return
1116
1117    # create a list of the phase variables
1118    Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtable,BLtable,MFtable,maxSSwave = G2stIO.GetPhaseData(Phases,rbIds=rbIds,Print=False)
1119    phaseList = []
1120    for item in phaseDict:
1121        if item.split(':')[2] not in badPhaseParms:
1122            phaseList.append(item)
1123    phaseList.sort()
1124    phaseAtNames = {}
1125    phaseAtTypes = {}
1126    TypeList = []
1127    for item in phaseList:
1128        Split = item.split(':')
1129        if Split[2][:2] in ['AU','Af','dA','AM']:
1130            Id = int(Split[0])
1131            phaseAtNames[item] = AtomDict[Id][int(Split[3])][0]
1132            phaseAtTypes[item] = AtomDict[Id][int(Split[3])][1]
1133            if phaseAtTypes[item] not in TypeList:
1134                TypeList.append(phaseAtTypes[item])
1135        else:
1136            phaseAtNames[item] = ''
1137            phaseAtTypes[item] = ''
1138             
1139    # create a list of the hist*phase variables
1140    seqList = G2frame.testSeqRefineMode()
1141    if seqList: # for sequential refinement, only process 1st histgram in list
1142        histDict = {seqList[0]:Histograms[seqList[0]]}
1143    else:
1144        histDict = Histograms
1145    hapVary,hapDict,controlDict = G2stIO.GetHistogramPhaseData(Phases,histDict,Print=False,resetRefList=False)
1146    hapList = sorted([i for i in hapDict.keys() if i.split(':')[2] not in ('Type',)])
1147    if seqList: # convert histogram # to wildcard
1148        wildList = [] # list of variables with "*" for histogram number
1149        for i in hapList:
1150            s = i.split(':')
1151            if s[1] == "": continue
1152            s[1] = '*'
1153            sj = ':'.join(s)
1154            if sj not in wildList: wildList.append(sj)
1155        hapList = wildList
1156    histVary,histDict,controlDict = G2stIO.GetHistogramData(histDict,Print=False)
1157    histList = list(histDict.keys())
1158    histList.sort()
1159    if seqList: # convert histogram # to wildcard
1160        wildList = [] # list of variables with "*" for histogram number
1161        for i in histList:
1162            s = i.split(':')
1163            if s[1] == "": continue
1164            s[1] = '*'
1165            sj = ':'.join(s)
1166            if sj not in wildList: wildList.append(sj)
1167        histList = wildList       
1168    Indx = {}
1169    G2frame.Page = [0,'phs']
1170   
1171    G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.ConstraintMenu)
1172    SetStatusLine('')
1173   
1174    G2frame.Bind(wx.EVT_MENU, OnAddConstraint, id=G2G.wxID_CONSTRAINTADD)
1175    G2frame.Bind(wx.EVT_MENU, OnAddFunction, id=G2G.wxID_FUNCTADD)
1176    G2frame.Bind(wx.EVT_MENU, OnAddEquivalence, id=G2G.wxID_EQUIVADD)
1177    G2frame.Bind(wx.EVT_MENU, OnAddHold, id=G2G.wxID_HOLDADD)
1178    G2frame.Bind(wx.EVT_MENU, OnAddAtomEquiv, id=G2G.wxID_EQUIVALANCEATOMS)
1179#    G2frame.Bind(wx.EVT_MENU, OnAddRiding, id=G2G.wxID_ADDRIDING)
1180    # tab commands
1181    for id in (G2G.wxID_CONSPHASE,
1182               G2G.wxID_CONSHAP,
1183               G2G.wxID_CONSHIST,
1184               G2G.wxID_CONSGLOBAL,
1185               G2G.wxID_CONSSYM,
1186               ):
1187        G2frame.Bind(wx.EVT_MENU, RaisePage,id=id)
1188
1189    #G2frame.constr = G2G.GSNoteBook(parent=G2frame.dataWindow,size=G2frame.dataWindow.GetClientSize())
1190    G2frame.constr = G2G.GSNoteBook(parent=G2frame.dataWindow)
1191    G2frame.dataWindow.GetSizer().Add(G2frame.constr,1,wx.ALL|wx.EXPAND)
1192    # note that order of pages is hard-coded in RaisePage
1193    PhaseConstr = wx.ScrolledWindow(G2frame.constr)
1194    G2frame.constr.AddPage(PhaseConstr,'Phase')
1195    HAPConstr = wx.ScrolledWindow(G2frame.constr)
1196    G2frame.constr.AddPage(HAPConstr,'Histogram/Phase')
1197    HistConstr = wx.ScrolledWindow(G2frame.constr)
1198    G2frame.constr.AddPage(HistConstr,'Histogram')
1199    GlobalConstr = wx.ScrolledWindow(G2frame.constr)
1200    G2frame.constr.AddPage(GlobalConstr,'Global')
1201    SymConstr = wx.ScrolledWindow(G2frame.constr)
1202    G2frame.constr.AddPage(SymConstr,'Sym-Generated')
1203    wx.CallAfter(OnPageChanged,None)
1204    G2frame.constr.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, OnPageChanged)
1205    # validate all the constrants -- should not see any errors here normally
1206    errmsg,warnmsg = ConstraintsCheck(data)
1207    if errmsg:
1208        G2frame.ErrorDialog('Constraint Error',
1209                            'Error in constraints:\n'+errmsg+'\nCheck console output for more information',
1210                            parent=G2frame)
1211        print (errmsg)
1212        print (G2mv.VarRemapShow([],True))
1213    elif warnmsg:
1214        print ('Unexpected contraint warning:\n'+warnmsg)
1215       
1216################################################################################
1217# check scale & phase fractions, create constraint if needed
1218################################################################################
1219def CheckAllScalePhaseFractions(G2frame):
1220    '''Check if scale factor and all phase fractions are refined without a constraint
1221    for all used histograms, if so, offer the user a chance to create a constraint
1222    on the sum of phase fractions
1223    '''
1224    histograms, phases = G2frame.GetUsedHistogramsAndPhasesfromTree()
1225    for i,hist in enumerate(histograms):
1226        CheckScalePhaseFractions(G2frame,hist,histograms,phases)
1227       
1228def CheckScalePhaseFractions(G2frame,hist,histograms,phases):
1229    '''Check if scale factor and all phase fractions are refined without a constraint
1230    for histogram hist, if so, offer the user a chance to create a constraint
1231    on the sum of phase fractions
1232    '''
1233    if G2frame.testSeqRefineMode():
1234        histStr = '*'
1235    else:
1236        histStr = str(histograms[hist]['hId'])
1237    # Is this powder?
1238    if not hist.startswith('PWDR '): return
1239    # do this only if the scale factor is varied
1240    if not histograms[hist]['Sample Parameters']['Scale'][1]: return
1241    # are all phase fractions varied in all used histograms?
1242    phaseCount = 0
1243    for p in phases:
1244        if hist not in phases[p]['Histograms']: continue
1245        if phases[p]['Histograms'][hist]['Use'] and not phases[p]['Histograms'][hist]['Scale'][1]:
1246            return
1247        else:
1248            phaseCount += 1
1249   
1250    # all phase fractions and scale factor varied, now scan through constraints
1251    sub = G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Constraints') 
1252    Constraints = G2frame.GPXtree.GetItemPyData(sub)
1253    for c in Constraints.get('HAP',[]):
1254        if c[-1] != 'c': continue
1255        if not c[-3]: continue
1256        if len(c[:-3]) != phaseCount: continue
1257        # got a constraint equation with right number of terms, is it on phase fractions for
1258        # the correct histogram?
1259        if all([(i[1].name == 'Scale' and i[1].varname().split(':')[1] == histStr) for i in c[:-3]]):
1260            # got a constraint, this is OK
1261            return
1262    dlg = wx.MessageDialog(G2frame,'You are refining the scale factor and all phase fractions for histogram #'+
1263        histStr+'. This will produce an unstable refinement. '+
1264        'Do you want to constrain the sum of phase fractions?','Create constraint?',wx.OK|wx.CANCEL)
1265    if dlg.ShowModal() != wx.ID_OK:
1266        dlg.Destroy()
1267        return
1268    dlg.Destroy()
1269
1270    constr = []
1271    for p in phases:
1272        if hist not in phases[p]['Histograms']: continue
1273        if not phases[p]['Histograms'][hist]['Use']: continue
1274        constr += [[1.0,G2obj.G2VarObj(':'.join((str(phases[p]['pId']),histStr,'Scale')))]]
1275    constr += [1.0,None,'c']
1276    Constraints['HAP'] += [constr]
1277       
1278################################################################################
1279#### Make nuclear/magnetic phase transition constraints - called by OnTransform in G2phsGUI
1280################################################################################       
1281       
1282def TransConstraints(G2frame,oldPhase,newPhase,Trans,Vec,atCodes):
1283    '''Add constraints for new magnetic phase created via transformation of old
1284    nuclear one
1285    NB: A = [G11,G22,G33,2*G12,2*G13,2*G23]
1286    '''
1287   
1288    def SetUniqAj(pId,iA,SGData):
1289        SGLaue = SGData['SGLaue']
1290        SGUniq = SGData['SGUniq']
1291        if SGLaue in ['m3','m3m']:
1292            if iA in [0,1,2]:
1293                parm = '%d::%s'%(pId,'A0')
1294            else:
1295                parm = None
1296        elif SGLaue in ['4/m','4/mmm']:
1297            if iA in [0,1]:
1298                parm = '%d::%s'%(pId,'A0')
1299            elif iA == 2:
1300                parm = '%d::%s'%(pId,'A2')
1301            else:
1302                parm = None
1303        elif SGLaue in ['6/m','6/mmm','3m1', '31m', '3']:
1304            if iA in [0,1,3]:
1305                parm = '%d::%s'%(pId,'A0')
1306            elif iA == 2:
1307                parm = '%d::%s'%(pId,'A2')
1308            else:
1309                parm = None
1310        elif SGLaue in ['3R', '3mR']:
1311            if ia in [0,1,2]:
1312                parm = '%d::%s'%(pId,'A0')
1313            else:
1314                parm = '%d::%s'%(pId,'A3')
1315        elif SGLaue in ['mmm',]:
1316            if iA in [0,1,2]:
1317                parm = '%d::A%s'%(pId,iA)
1318            else:
1319                parm = None
1320        elif SGLaue == '2/m':
1321            if iA in [0,1,2]:
1322                parm = '%d::A%s'%(pId,iA)
1323            elif iA == 3 and SGUniq == 'c':
1324                parm = '%d::A%s'%(pId,iA)
1325            elif iA == 4 and SGUniq == 'b':
1326                parm = '%d::A%s'%(pId,iA)
1327            elif iA == 5 and SGUniq == 'a':
1328                parm = '%d::A%s'%(pId,iA)
1329            else:
1330                parm = None           
1331        else:
1332            parm = '%d::A%s'%(pId,iA)
1333        return parm
1334   
1335    Histograms,Phases = G2frame.GetUsedHistogramsAndPhasesfromTree()
1336    UseList = newPhase['Histograms']
1337    detTrans = np.abs(nl.det(Trans))
1338    opId = oldPhase['pId']
1339    npId = newPhase['pId']
1340    cx,ct,cs,cia = newPhase['General']['AtomPtrs']
1341    nAtoms = newPhase['Atoms']
1342    oSGData = oldPhase['General']['SGData']
1343    nSGData = newPhase['General']['SGData']
1344    #oAcof = G2lat.cell2A(oldPhase['General']['Cell'][1:7])
1345    #nAcof = G2lat.cell2A(newPhase['General']['Cell'][1:7])
1346    item = G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Constraints')
1347    if not item:
1348        print('Error: no constraints in Data Tree')
1349        return
1350    constraints = G2frame.GPXtree.GetItemPyData(item)
1351    xnames = ['dAx','dAy','dAz']
1352    # constraints on matching atom params between phases
1353    for ia,code in enumerate(atCodes):
1354        atom = nAtoms[ia]
1355        if not ia and atom[cia] == 'A':
1356            wx.MessageDialog(G2frame,
1357                'Anisotropic thermal motion constraints are not developed at the present time',
1358                'Anisotropic thermal constraint?',style=wx.ICON_INFORMATION).ShowModal()
1359        siteSym = G2spc.SytSym(atom[cx:cx+3],nSGData)[0]
1360        CSX = G2spc.GetCSxinel(siteSym)
1361#        CSU = G2spc.GetCSuinel(siteSym)
1362        item = code.split('+')[0]
1363        iat,opr = item.split(':')
1364        Nop = abs(int(opr))%100-1
1365        if '-' in opr:
1366            Nop *= -1
1367        Opr = oldPhase['General']['SGData']['SGOps'][abs(Nop)][0]
1368        if Nop < 0:         #inversion
1369            Opr *= -1
1370        XOpr = np.inner(Opr,Trans)
1371        for ix in list(set(CSX[0])):
1372            if not ix:
1373                continue
1374            name = xnames[ix-1]
1375            IndpCon = [1.0,G2obj.G2VarObj('%d::%s:%d'%(npId,name,ia))]
1376            DepCons = []
1377            for iop,opval in enumerate(XOpr[ix-1]):
1378                if opval:
1379                    DepCons.append([opval,G2obj.G2VarObj('%d::%s:%s'%(opId,xnames[iop],iat))])
1380            if len(DepCons) == 1:
1381                constraints['Phase'].append([DepCons[0],IndpCon,None,None,'e'])
1382            elif len(DepCons) > 1:
1383                for Dep in DepCons:
1384                    Dep[0] *= -1
1385                constraints['Phase'].append([IndpCon]+DepCons+[0.0,None,'c'])
1386        for name in ['Afrac','AUiso']:
1387            IndpCon = [1.0,G2obj.G2VarObj('%d::%s:%d'%(npId,name,ia))]
1388            DepCons = [1.0,G2obj.G2VarObj('%d::%s:%s'%(opId,name,iat))]
1389            constraints['Phase'].append([DepCons,IndpCon,None,None,'e'])
1390           
1391        # unfinished Anisotropic constraint generation
1392#        Uids = [[0,0,'AU11'],[1,1,'AU22'],[2,2,'AU33'],[0,1,'AU12'],[0,2,'AU13'],[1,2,'AU23']]
1393#        DepConsDict = dict(zip(Us,[[],[],[],[],[],[]]))
1394#        for iu,Uid in enumerate(Uids):
1395#            UMT = np.zeros((3,3))
1396#            UMT[Uid[0],Uid[1]] = 1
1397#            nUMT = G2lat.prodMGMT(UMT,invTrans)
1398#            nUT = G2lat.UijtoU6(nUMT)
1399#            for iu,nU in enumerate(nUT):
1400#                if abs(nU) > 1.e-8:
1401#                    parm = '%d::%s;%s'%(opId,Us[iu],iat)
1402#                    DepConsDict[Uid[2]].append([abs(nU%1.),G2obj.G2VarObj(parm)])
1403#        nUcof = atom[iu:iu+6]
1404#        conStrings = []
1405#        for iU,Usi in enumerate(Us):
1406#            parm = '%d::%s;%d'%(npId,Usi,ia)
1407#            parmDict[parm] = nUcof[iU]
1408#            varyList.append(parm)
1409#            IndpCon = [1.0,G2obj.G2VarObj(parm)]
1410#            conStr = str([IndpCon,DepConsDict[Usi]])
1411#            if conStr in conStrings:
1412#                continue
1413#            conStrings.append(conStr)
1414#            if len(DepConsDict[Usi]) == 1:
1415#                if DepConsDict[Usi][0]:
1416#                    constraints['Phase'].append([IndpCon,DepConsDict[Usi][0],None,None,'e'])
1417#            elif len(DepConsDict[Usi]) > 1:       
1418#                for Dep in DepConsDict[Usi]:
1419#                    Dep[0] *= -1
1420#                constraints['Phase'].append([IndpCon]+DepConsDict[Usi]+[0.0,None,'c'])
1421           
1422        #how do I do Uij's for most Trans?
1423
1424    # constraints on lattice parameters between phases
1425    T = nl.inv(Trans).T
1426    conMat = [
1427        [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]],
1428        [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]],
1429        [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]],
1430        [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]],
1431        [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]],
1432        [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]]
1433        ]
1434    # Gnew = conMat * A:
1435#         T00**2*a0  T01**2*a1 T02**2*a2 T00*T01*a3    T00*T02*a4    T01*T02*a5
1436#         T10**2*a0  T11**2*a1 T12**2*a2 T10*T11*a3    T10*T12*a4    T11*T12*a5
1437#         T20**2*a0  T21**2*a1 T22**2*a2 T20*T21*a3    T20*T22*a4    T21*T22*a5
1438#         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
1439#         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
1440#         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
1441    # Generated as symbolic code using:
1442    # import sym
1443    # A0, A1, A2, A3, A4, A5 = sym.symbols('A0, A1, A2, A3, A4, A5')
1444    # G = sym.Matrix([ [A0,    A3/2.,  A4/2.], [A3/2.,  A1,    A5/2.], [A4/2., A5/2.,    A2]])
1445    # transformation matrix
1446    # T00, T10, T20, T01, T11, T21, T02, T12, T22 = sym.symbols('T00, T10, T20, T01, T11, T21, T02, T12, T22')
1447    # Tr = sym.Matrix([ [T00, T10, T20], [T01, T11, T21], [T02, T12, T22],])
1448    # Gnew = (Tr.T*G)*Tr
1449   
1450    #print('old A',G2lat.cell2A(oldPhase['General']['Cell'][1:7]))
1451    #print('new A',G2lat.cell2A(newPhase['General']['Cell'][1:7]))
1452    for iAnew,Asi in enumerate(['A0','A1','A2','A3','A4','A5']): # loop through A[i] for new cell
1453        Nparm = str(npId) + '::' + Asi
1454        if Nparm != SetUniqAj(npId,iAnew,nSGData):
1455            continue # skip: Ai constrained from Aj or must be zero
1456        multDict = {}
1457        for iAorg in range(6):
1458            cA = conMat[iAnew][iAorg] # coeff for A[i] in constraint matrix
1459            if abs(cA) < 1.e-8: continue
1460            parm = SetUniqAj(opId,iAorg,oSGData) # translate to unique A[i] in original cell
1461            if not parm: continue # must be zero
1462            # sum coeff
1463            if parm in multDict:
1464                multDict[parm] += cA
1465            else:
1466                multDict[parm] = cA
1467        # any non-zero multipliers?
1468        maxMult = 0
1469        for i in multDict:
1470            maxMult = max(maxMult,abs(multDict[i]))
1471        if maxMult <= 0:  # Nparm computes as zero; Fix this parameter
1472            constraints['Phase'] += [[
1473                [0.0,G2obj.G2VarObj(Nparm)],
1474                None,None,'h']]
1475        elif len(multDict) == 1:        # create equivalence
1476            key = list(multDict.keys())[0]
1477            constraints['Phase'] += [[
1478                [1.0,G2obj.G2VarObj(key)],
1479                [multDict[key],G2obj.G2VarObj(Nparm)],
1480                None,None,'e']]
1481        else:                           # create constraint
1482            constr = [[-1.0,G2obj.G2VarObj(Nparm)]]
1483            for key in multDict:
1484                constr += [[multDict[key],G2obj.G2VarObj(key)]]
1485            constr += [0.0,None,'c']
1486            constraints['Phase'] += [constr]
1487   
1488    # constraints on HAP Scale, etc.
1489    for hId,hist in enumerate(UseList):    #HAP - seems OK
1490        ohapkey = '%d:%d:'%(opId,hId)
1491        nhapkey = '%d:%d:'%(npId,hId)
1492        IndpCon = [1.0,G2obj.G2VarObj(ohapkey+'Scale')]
1493        DepCons = [detTrans,G2obj.G2VarObj(nhapkey+'Scale')]
1494        constraints['HAP'].append([DepCons,IndpCon,None,None,'e'])
1495        for name in ['Size;i','Mustrain;i']:
1496            IndpCon = [1.0,G2obj.G2VarObj(ohapkey+name)]
1497            DepCons = [1.0,G2obj.G2VarObj(nhapkey+name)]
1498            constraints['HAP'].append([IndpCon,DepCons,None,None,'e'])
1499       
1500################################################################################
1501#### Rigid bodies
1502################################################################################
1503
1504def UpdateRigidBodies(G2frame,data):
1505    '''Called when Rigid bodies tree item is selected.
1506    Displays the rigid bodies in the data window
1507    '''
1508    if not data.get('RBIds') or not data:
1509        data.update({'Vector':{'AtInfo':{}},'Residue':{'AtInfo':{}},
1510            'RBIds':{'Vector':[],'Residue':[]}})       #empty/bad dict - fill it
1511           
1512    global resList,rbId
1513    Indx = {}
1514    resList = []
1515    plotDefaults = {'oldxy':[0.,0.],'Quaternion':[0.,0.,0.,1.],'cameraPos':30.,'viewDir':[0,0,1],}
1516
1517    G2frame.rbBook = G2G.GSNoteBook(parent=G2frame.dataWindow)
1518    G2frame.dataWindow.GetSizer().Add(G2frame.rbBook,1,wx.ALL|wx.EXPAND)
1519    VectorRB = wx.ScrolledWindow(G2frame.rbBook)
1520    VectorRBDisplay = wx.Panel(VectorRB)
1521    G2frame.rbBook.AddPage(VectorRB,'Vector rigid bodies')
1522    ResidueRB = wx.ScrolledWindow(G2frame.rbBook)
1523    ResidueRBDisplay = wx.Panel(ResidueRB)
1524    G2frame.rbBook.AddPage(ResidueRB,'Residue rigid bodies')
1525   
1526    def OnPageChanged(event):
1527        global resList
1528        resList = []
1529        if event:       #page change event!
1530            page = event.GetSelection()
1531        else:
1532            page = G2frame.rbBook.GetSelection()
1533        G2frame.rbBook.ChangeSelection(page)
1534        text = G2frame.rbBook.GetPageText(page)
1535        if text == 'Vector rigid bodies':
1536            G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.VectorBodyMenu)
1537            G2frame.Bind(wx.EVT_MENU, AddVectorRB, id=G2G.wxID_VECTORBODYADD)
1538            G2frame.Page = [page,'vrb']
1539            UpdateVectorRB()
1540        elif text == 'Residue rigid bodies':
1541            G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.RigidBodyMenu)
1542            G2frame.Bind(wx.EVT_MENU, AddResidueRB, id=G2G.wxID_RIGIDBODYADD)
1543            G2frame.Bind(wx.EVT_MENU, OnImportRigidBody, id=G2G.wxID_RIGIDBODYIMPORT)
1544            G2frame.Bind(wx.EVT_MENU, OnDefineTorsSeq, id=G2G.wxID_RESIDUETORSSEQ) #enable only if residue RBs exist?
1545            G2frame.Page = [page,'rrb']
1546            UpdateResidueRB()
1547           
1548    def getMacroFile(macName):
1549        defDir = os.path.join(os.path.split(__file__)[0],'GSASIImacros')
1550        dlg = wx.FileDialog(G2frame,message='Choose '+macName+' rigid body macro file',
1551            defaultDir=defDir,defaultFile="",wildcard="GSAS-II macro file (*.mac)|*.mac",
1552            style=wx.FD_OPEN | wx.FD_CHANGE_DIR)
1553        try:
1554            if dlg.ShowModal() == wx.ID_OK:
1555                macfile = dlg.GetPath()
1556                macro = open(macfile,'Ur')
1557                head = macro.readline()
1558                if macName not in head:
1559                    print (head)
1560                    print ('**** ERROR - wrong restraint macro file selected, try again ****')
1561                    macro = []
1562            else: # cancel was pressed
1563                macro = []
1564        finally:
1565            dlg.Destroy()
1566        return macro        #advanced past 1st line
1567       
1568    def getTextFile():
1569        dlg = wx.FileDialog(G2frame,'Choose rigid body text file', '.', '',
1570            "GSAS-II text file (*.txt)|*.txt|XYZ file (*.xyz)|*.xyz|"
1571            "Sybyl mol2 file (*.mol2)|*.mol2|PDB file (*.pdb;*.ent)|*.pdb;*.ent",
1572            wx.FD_OPEN | wx.FD_CHANGE_DIR)
1573        try:
1574            if dlg.ShowModal() == wx.ID_OK:
1575                txtfile = dlg.GetPath()
1576                ext = os.path.splitext(txtfile)[1]
1577                text = open(txtfile,'Ur')
1578            else: # cancel was pressed
1579                ext = ''
1580                text = []
1581        finally:
1582            dlg.Destroy()
1583        if 'ent' in ext:
1584            ext = '.pdb'
1585        return text,ext.lower()
1586       
1587    def OnImportRigidBody(event):
1588        page = G2frame.rbBook.GetSelection()
1589        if 'Vector' in G2frame.rbBook.GetPageText(page):
1590            pass
1591        elif 'Residue' in G2frame.rbBook.GetPageText(page):
1592            ImportResidueRB()
1593           
1594    def AddVectorRB(event):
1595        AtInfo = data['Vector']['AtInfo']
1596        dlg = G2G.MultiIntegerDialog(G2frame,'New Rigid Body',['No. atoms','No. translations'],[1,1])
1597        if dlg.ShowModal() == wx.ID_OK:
1598            nAtoms,nTrans = dlg.GetValues()
1599            rbId = ran.randint(0,sys.maxsize)
1600            vecMag = [1.0 for i in range(nTrans)]
1601            vecRef = [False for i in range(nTrans)]
1602            vecVal = [np.zeros((nAtoms,3)) for j in range(nTrans)]
1603            rbTypes = ['C' for i in range(nAtoms)]
1604            Info = G2elem.GetAtomInfo('C')
1605            AtInfo['C'] = [Info['Drad'],Info['Color']]
1606            data['Vector'][rbId] = {'RBname':'UNKRB','VectMag':vecMag,'rbXYZ':np.zeros((nAtoms,3)),
1607                'rbRef':[0,1,2,False],'VectRef':vecRef,'rbTypes':rbTypes,'rbVect':vecVal,'useCount':0}
1608            data['RBIds']['Vector'].append(rbId)
1609        dlg.Destroy()
1610        UpdateVectorRB()
1611       
1612    def AddResidueRB(event):
1613        global rbId
1614        AtInfo = data['Residue']['AtInfo']
1615        macro = getMacroFile('rigid body')
1616        if not macro:
1617            return
1618        macStr = macro.readline()
1619        while macStr:
1620            items = macStr.split()
1621            if 'I' == items[0]:
1622                rbId = ran.randint(0,sys.maxsize)
1623                rbName = items[1]
1624                rbTypes = []
1625                rbXYZ = []
1626                rbSeq = []
1627                atNames = []
1628                nAtms,nSeq,nOrig,mRef,nRef = [int(items[i]) for i in [2,3,4,5,6]]
1629                for iAtm in range(nAtms):
1630                    macStr = macro.readline().split()
1631                    atName = macStr[0]
1632                    atType = macStr[1]
1633                    atNames.append(atName)
1634                    rbXYZ.append([float(macStr[i]) for i in [2,3,4]])
1635                    rbTypes.append(atType)
1636                    if atType not in AtInfo:
1637                        Info = G2elem.GetAtomInfo(atType)
1638                        AtInfo[atType] = [Info['Drad'],Info['Color']]
1639                rbXYZ = np.array(rbXYZ)-np.array(rbXYZ[nOrig-1])
1640                for iSeq in range(nSeq):
1641                    macStr = macro.readline().split()
1642                    mSeq = int(macStr[0])
1643                    for jSeq in range(mSeq):
1644                        macStr = macro.readline().split()
1645                        iBeg = int(macStr[0])-1
1646                        iFin = int(macStr[1])-1
1647                        angle = 0.0
1648                        nMove = int(macStr[2])
1649                        iMove = [int(macStr[i])-1 for i in range(3,nMove+3)]
1650                        rbSeq.append([iBeg,iFin,angle,iMove])
1651                data['Residue'][rbId] = {'RBname':rbName,'rbXYZ':rbXYZ,'rbTypes':rbTypes,
1652                    'atNames':atNames,'rbRef':[nOrig-1,mRef-1,nRef-1,True],'rbSeq':rbSeq,
1653                    'SelSeq':[0,0],'useCount':0}
1654                data['RBIds']['Residue'].append(rbId)
1655                print ('Rigid body '+rbName+' added')
1656            macStr = macro.readline()
1657        macro.close()
1658        UpdateResidueRB()
1659       
1660    def ImportResidueRB():
1661        global rbId
1662        AtInfo = data['Residue']['AtInfo']
1663        text,ext = getTextFile()
1664        if not text:
1665            return
1666        rbId = ran.randint(0,sys.maxsize)
1667        rbTypes = []
1668        rbXYZ = []
1669        atNames = []
1670        txtStr = text.readline()
1671        if 'xyz' in ext:
1672            txtStr = text.readline()
1673            txtStr = text.readline()
1674        elif 'mol2' in ext:
1675            while 'ATOM' not in txtStr:
1676                txtStr = text.readline()
1677            txtStr = text.readline()
1678        elif 'pdb' in ext:
1679            while 'ATOM' not in txtStr[:6] and 'HETATM' not in txtStr[:6]:
1680                txtStr = text.readline()
1681                #print txtStr
1682        items = txtStr.split()
1683        while len(items):
1684            if 'txt' in ext:
1685                atName = items[0]
1686                atType = items[1]
1687                rbXYZ.append([float(items[i]) for i in [2,3,4]])
1688            elif 'xyz' in ext:
1689                atType = items[0]
1690                rbXYZ.append([float(items[i]) for i in [1,2,3]])
1691                atName = atType+str(len(rbXYZ))
1692            elif 'mol2' in ext:
1693                atType = items[1]
1694                atName = items[1]+items[0]
1695                rbXYZ.append([float(items[i]) for i in [2,3,4]])
1696            elif 'pdb' in ext:
1697                atType = items[-1]
1698                atName = items[2]
1699                xyz = txtStr[30:55].split()                   
1700                rbXYZ.append([float(x) for x in xyz])
1701            atNames.append(atName)
1702            rbTypes.append(atType)
1703            if atType not in AtInfo:
1704                Info = G2elem.GetAtomInfo(atType)
1705                AtInfo[atType] = [Info['Drad'],Info['Color']]
1706            txtStr = text.readline()
1707            if 'mol2' in ext and 'BOND' in txtStr:
1708                break
1709            if 'pdb' in ext and ('ATOM' not in txtStr[:6] and 'HETATM' not in txtStr[:6]):
1710                break
1711            items = txtStr.split()
1712        if len(atNames) < 3:
1713            G2G.G2MessageBox(G2frame,'Not enough atoms in rigid body; must be 3 or more')
1714        else:
1715            rbXYZ = np.array(rbXYZ)-np.array(rbXYZ[0])
1716            Xxyz = rbXYZ[1]
1717            X = Xxyz/np.sqrt(np.sum(Xxyz**2))
1718            Yxyz = rbXYZ[2]
1719            Y = Yxyz/np.sqrt(np.sum(Yxyz**2))
1720            Mat = G2mth.getRBTransMat(X,Y)
1721            rbXYZ = np.inner(Mat,rbXYZ).T
1722            data['Residue'][rbId] = {'RBname':'UNKRB','rbXYZ':rbXYZ,'rbTypes':rbTypes,
1723                'atNames':atNames,'rbRef':[0,1,2,False],'rbSeq':[],'SelSeq':[0,0],'useCount':0}
1724            data['RBIds']['Residue'].append(rbId)
1725            print ('Rigid body UNKRB added')
1726        text.close()
1727        UpdateResidueRB(rbId)
1728       
1729    def FindNeighbors(Orig,XYZ,atTypes,atNames,AtInfo):
1730        Radii = []
1731        for Atype in atTypes:
1732            Radii.append(AtInfo[Atype][0])
1733        Radii = np.array(Radii)
1734        Neigh = []
1735        Dx = XYZ-XYZ[Orig]
1736        dist = np.sqrt(np.sum(Dx**2,axis=1))
1737        sumR = Radii[Orig]+Radii
1738        IndB = ma.nonzero(ma.masked_greater(dist-0.85*sumR,0.))
1739        for j in IndB[0]:
1740            if j != Orig and atTypes[j] != 'H':
1741                Neigh.append(atNames[j])
1742        return Neigh
1743       
1744    def FindAllNeighbors(XYZ,atTypes,atNames,AtInfo):
1745        NeighDict = {}
1746        for iat,xyz in enumerate(atNames):
1747            NeighDict[atNames[iat]] = FindNeighbors(iat,XYZ,atTypes,atNames,AtInfo)
1748        return NeighDict
1749       
1750    def FindRiding(Orig,Pivot,NeighDict):
1751        riding = [Orig,Pivot]
1752        iAdd = 1
1753        new = True
1754        while new:
1755            newAtms = NeighDict[riding[iAdd]]
1756            for At in newAtms:
1757                new = False
1758                if At not in riding:
1759                    riding.append(At)
1760                    new = True
1761            iAdd += 1
1762            if iAdd < len(riding):
1763                new = True
1764        return riding[2:]
1765                       
1766    def OnDefineTorsSeq(event):
1767        global rbId
1768        rbData = data['Residue'][rbId]
1769        if not len(rbData):
1770            return
1771        atNames = rbData['atNames']
1772        AtInfo = data['Residue']['AtInfo']
1773        atTypes = rbData['rbTypes']
1774        XYZ = rbData['rbXYZ']
1775        neighDict = FindAllNeighbors(XYZ,atTypes,atNames,AtInfo)
1776        TargList = []           
1777        dlg = wx.SingleChoiceDialog(G2frame,'Select origin atom for torsion sequence','Origin atom',rbData['atNames'])
1778        if dlg.ShowModal() == wx.ID_OK:
1779            Orig = dlg.GetSelection()
1780            TargList = neighDict[atNames[Orig]]
1781        dlg.Destroy()
1782        if not len(TargList):
1783            return
1784        dlg = wx.SingleChoiceDialog(G2frame,'Select pivot atom for torsion sequence','Pivot atom',TargList)
1785        if dlg.ShowModal() == wx.ID_OK:
1786            Piv = atNames.index(TargList[dlg.GetSelection()])
1787            riding = FindRiding(atNames[Orig],atNames[Piv],neighDict)
1788            Riding = []
1789            for atm in riding:
1790                Riding.append(atNames.index(atm))
1791            rbData['rbSeq'].append([Orig,Piv,0.0,Riding])           
1792        dlg.Destroy()
1793        UpdateResidueRB(rbId)
1794
1795    def UpdateVectorRB(Scroll=0):
1796        AtInfo = data['Vector']['AtInfo']
1797        refChoice = {}
1798        if 'DELETED' in str(G2frame.GetStatusBar()):   #seems to be no other way to do this (wx bug)
1799            if GSASIIpath.GetConfigValue('debug'):
1800                print ('DBG_wx error: Rigid Body/Status not cleanly deleted after Refine')
1801            return
1802        SetStatusLine(' You may use e.g. "c60" or "s60" for a vector entry')
1803        def rbNameSizer(rbId,rbData):
1804
1805            def OnRBName(event):
1806                event.Skip()
1807                Obj = event.GetEventObject()
1808                rbData['RBname'] = Obj.GetValue()
1809               
1810            def OnDelRB(event):
1811                Obj = event.GetEventObject()
1812                rbId = Indx[Obj.GetId()]
1813                if rbId in data['Vector']:
1814                    del data['Vector'][rbId]
1815                    data['RBIds']['Vector'].remove(rbId)
1816                    rbData['useCount'] -= 1
1817                wx.CallAfter(UpdateVectorRB)
1818               
1819            def OnPlotRB(event):
1820                Obj = event.GetEventObject()
1821                Obj.SetValue(False)
1822                G2plt.PlotRigidBody(G2frame,'Vector',AtInfo,rbData,plotDefaults)
1823           
1824            nameSizer = wx.BoxSizer(wx.HORIZONTAL)
1825            nameSizer.Add(wx.StaticText(VectorRBDisplay,-1,'Rigid body name: '),
1826                0,wx.ALIGN_CENTER_VERTICAL)
1827            RBname = wx.TextCtrl(VectorRBDisplay,-1,rbData['RBname'])
1828            Indx[RBname.GetId()] = rbId
1829            RBname.Bind(wx.EVT_TEXT_ENTER,OnRBName)
1830            RBname.Bind(wx.EVT_KILL_FOCUS,OnRBName)
1831            nameSizer.Add(RBname,0,wx.ALIGN_CENTER_VERTICAL)
1832            nameSizer.Add((5,0),)
1833            plotRB = wx.CheckBox(VectorRBDisplay,-1,'Plot?')
1834            Indx[plotRB.GetId()] = rbId
1835            plotRB.Bind(wx.EVT_CHECKBOX,OnPlotRB)
1836            nameSizer.Add(plotRB,0,wx.ALIGN_CENTER_VERTICAL)
1837            nameSizer.Add((5,0),)
1838            if not rbData['useCount']:
1839                delRB = wx.CheckBox(VectorRBDisplay,-1,'Delete?')
1840                Indx[delRB.GetId()] = rbId
1841                delRB.Bind(wx.EVT_CHECKBOX,OnDelRB)
1842                nameSizer.Add(delRB,0,wx.ALIGN_CENTER_VERTICAL)
1843            return nameSizer
1844           
1845        def rbRefAtmSizer(rbId,rbData):
1846           
1847            def OnRefSel(event):
1848                Obj = event.GetEventObject()
1849                iref = Indx[Obj.GetId()]
1850                sel = Obj.GetValue()
1851                rbData['rbRef'][iref] = atNames.index(sel)
1852                FillRefChoice(rbId,rbData)
1853           
1854            refAtmSizer = wx.BoxSizer(wx.HORIZONTAL)
1855            atNames = [name+str(i) for i,name in enumerate(rbData['rbTypes'])]
1856            rbRef = rbData.get('rbRef',[0,1,2,False])
1857            rbData['rbRef'] = rbRef
1858            if rbData['useCount']:
1859                refAtmSizer.Add(wx.StaticText(VectorRBDisplay,-1,
1860                    'Orientation reference atoms A-B-C: %s, %s, %s'%(atNames[rbRef[0]], \
1861                     atNames[rbRef[1]],atNames[rbRef[2]])),0)
1862            else:
1863                refAtmSizer.Add(wx.StaticText(VectorRBDisplay,-1,
1864                    'Orientation reference atoms A-B-C: '),0,wx.ALIGN_CENTER_VERTICAL)
1865                for i in range(3):
1866                    choices = [atNames[j] for j in refChoice[rbId][i]]
1867                    refSel = wx.ComboBox(VectorRBDisplay,-1,value='',
1868                        choices=choices,style=wx.CB_READONLY|wx.CB_DROPDOWN)
1869                    refSel.SetValue(atNames[rbRef[i]])
1870                    refSel.Bind(wx.EVT_COMBOBOX, OnRefSel)
1871                    Indx[refSel.GetId()] = i
1872                    refAtmSizer.Add(refSel,0,wx.ALIGN_CENTER_VERTICAL)
1873            return refAtmSizer
1874                       
1875        def rbVectMag(rbId,imag,rbData):
1876           
1877            def OnRBVectorMag(event):
1878                event.Skip()
1879                Obj = event.GetEventObject()
1880                rbId,imag = Indx[Obj.GetId()]
1881                try:
1882                    val = float(Obj.GetValue())
1883                    if val <= 0.:
1884                        raise ValueError
1885                    rbData['VectMag'][imag] = val
1886                except ValueError:
1887                    pass
1888                Obj.SetValue('%8.4f'%(val))
1889                wx.CallAfter(UpdateVectorRB,VectorRB.GetScrollPos(wx.VERTICAL))
1890                G2plt.PlotRigidBody(G2frame,'Vector',AtInfo,data['Vector'][rbId],plotDefaults)
1891               
1892            def OnRBVectorRef(event):
1893                Obj = event.GetEventObject()
1894                rbId,imag = Indx[Obj.GetId()]
1895                rbData['VectRef'][imag] = Obj.GetValue()
1896                       
1897            magSizer = wx.BoxSizer(wx.HORIZONTAL)
1898            magSizer.Add(wx.StaticText(VectorRBDisplay,-1,'Translation magnitude: '),
1899                0,wx.ALIGN_CENTER_VERTICAL)
1900            magValue = wx.TextCtrl(VectorRBDisplay,-1,'%8.4f'%(rbData['VectMag'][imag]))
1901            Indx[magValue.GetId()] = [rbId,imag]
1902            magValue.Bind(wx.EVT_TEXT_ENTER,OnRBVectorMag)
1903            magValue.Bind(wx.EVT_KILL_FOCUS,OnRBVectorMag)
1904            magSizer.Add(magValue,0,wx.ALIGN_CENTER_VERTICAL)
1905            magSizer.Add((5,0),)
1906            magref = wx.CheckBox(VectorRBDisplay,-1,label=' Refine?') 
1907            magref.SetValue(rbData['VectRef'][imag])
1908            magref.Bind(wx.EVT_CHECKBOX,OnRBVectorRef)
1909            Indx[magref.GetId()] = [rbId,imag]
1910            magSizer.Add(magref,0,wx.ALIGN_CENTER_VERTICAL)
1911            return magSizer
1912           
1913        def rbVectors(rbId,imag,mag,XYZ,rbData):
1914
1915            def TypeSelect(event):
1916                AtInfo = data['Vector']['AtInfo']
1917                r,c = event.GetRow(),event.GetCol()
1918                if vecGrid.GetColLabelValue(c) == 'Type':
1919                    PE = G2elemGUI.PickElement(G2frame,oneOnly=True)
1920                    if PE.ShowModal() == wx.ID_OK:
1921                        if PE.Elem != 'None':
1922                            El = PE.Elem.strip().lower().capitalize()
1923                            if El not in AtInfo:
1924                                Info = G2elem.GetAtomInfo(El)
1925                                AtInfo[El] = [Info['Drad'],Info['Color']]
1926                            rbData['rbTypes'][r] = El
1927                            vecGrid.SetCellValue(r,c,El)
1928                    PE.Destroy()
1929                wx.CallAfter(UpdateVectorRB,VectorRB.GetScrollPos(wx.VERTICAL))
1930
1931            def ChangeCell(event):
1932                r,c =  event.GetRow(),event.GetCol()
1933                if r >= 0 and (0 <= c < 3):
1934                    try:
1935                        val = float(vecGrid.GetCellValue(r,c))
1936                        rbData['rbVect'][imag][r][c] = val
1937                    except ValueError:
1938                        pass
1939                G2plt.PlotRigidBody(G2frame,'Vector',AtInfo,data['Vector'][rbId],plotDefaults)
1940                wx.CallAfter(UpdateVectorRB,VectorRB.GetScrollPos(wx.VERTICAL))
1941
1942            vecSizer = wx.BoxSizer()
1943            Types = 3*[wg.GRID_VALUE_FLOAT+':10,5',]+[wg.GRID_VALUE_STRING,]+3*[wg.GRID_VALUE_FLOAT+':10,5',]
1944            colLabels = ['Vector x','Vector y','Vector z','Type','Cart x','Cart y','Cart z']
1945            table = []
1946            rowLabels = []
1947            for ivec,xyz in enumerate(rbData['rbVect'][imag]):
1948                table.append(list(xyz)+[rbData['rbTypes'][ivec],]+list(XYZ[ivec]))
1949                rowLabels.append(str(ivec))
1950            vecTable = G2G.Table(table,rowLabels=rowLabels,colLabels=colLabels,types=Types)
1951            vecGrid = G2G.GSGrid(VectorRBDisplay)
1952            vecGrid.SetTable(vecTable, True)
1953            if 'phoenix' in wx.version():
1954                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGED, ChangeCell)
1955            else:
1956                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGE, ChangeCell)
1957            if not imag:
1958                vecGrid.Bind(wg.EVT_GRID_CELL_LEFT_DCLICK, TypeSelect)
1959            attr = wx.grid.GridCellAttr()
1960            attr.SetEditor(G2G.GridFractionEditor(vecGrid))
1961            for c in range(3):
1962                vecGrid.SetColAttr(c, attr)
1963            for row in range(vecTable.GetNumberRows()):
1964                if imag:
1965                    vecGrid.SetCellStyle(row,3,VERY_LIGHT_GREY,True)                   
1966                for col in [4,5,6]:
1967                    vecGrid.SetCellStyle(row,col,VERY_LIGHT_GREY,True)
1968#            vecGrid.SetScrollRate(0,0)
1969            vecGrid.AutoSizeColumns(False)
1970            vecSizer.Add(vecGrid)
1971            return vecSizer
1972       
1973        def FillRefChoice(rbId,rbData):
1974            choiceIds = [i for i in range(len(rbData['rbTypes']))]
1975           
1976            rbRef = rbData.get('rbRef',[-1,-1,-1,False])
1977            for i in range(3):
1978                choiceIds.remove(rbRef[i])
1979            refChoice[rbId] = [choiceIds[:],choiceIds[:],choiceIds[:]]
1980            for i in range(3):
1981                refChoice[rbId][i].append(rbRef[i])
1982                refChoice[rbId][i].sort()     
1983           
1984        if VectorRB.GetSizer(): VectorRB.GetSizer().Clear(True)
1985        VectorRBSizer = wx.BoxSizer(wx.VERTICAL)
1986        for rbId in data['RBIds']['Vector']:
1987            if rbId != 'AtInfo':
1988                rbData = data['Vector'][rbId]
1989                FillRefChoice(rbId,rbData)
1990                VectorRBSizer.Add(rbNameSizer(rbId,rbData),0)
1991                VectorRBSizer.Add(rbRefAtmSizer(rbId,rbData),0)
1992                XYZ = np.array([[0.,0.,0.] for Ty in rbData['rbTypes']])
1993                for imag,mag in enumerate(rbData['VectMag']):
1994                    XYZ += mag*rbData['rbVect'][imag]
1995                    VectorRBSizer.Add(rbVectMag(rbId,imag,rbData),0)
1996                    VectorRBSizer.Add(rbVectors(rbId,imag,mag,XYZ,rbData),0)
1997                VectorRBSizer.Add((5,5),0)
1998                data['Vector'][rbId]['rbXYZ'] = XYZ       
1999        VectorRBSizer.Layout()   
2000        VectorRBDisplay.SetSizer(VectorRBSizer,True)
2001        Size = VectorRBSizer.GetMinSize()
2002        Size[0] += 40
2003        Size[1] = max(Size[1],450) + 20
2004        VectorRBDisplay.SetSize(Size)
2005        VectorRB.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
2006        VectorRB.Scroll(0,Scroll)
2007       
2008    def UpdateResidueRB(rbId=0):
2009        AtInfo = data['Residue']['AtInfo']
2010        refChoice = {}
2011        RefObjs = []
2012
2013        def rbNameSizer(rbId,rbData):
2014
2015            def OnRBName(event):
2016                Obj = event.GetEventObject()
2017                rbData['RBname'] = Obj.GetValue()
2018                wx.CallAfter(UpdateResidueRB,rbId)
2019               
2020            def OnDelRB(event):
2021                Obj = event.GetEventObject()
2022                rbId = Indx[Obj.GetId()]
2023                if rbId in data['Residue']: 
2024                    del data['Residue'][rbId]
2025                    data['RBIds']['Residue'].remove(rbId)
2026                wx.CallAfter(UpdateResidueRB)
2027               
2028            def OnStripH(event):
2029                Obj = event.GetEventObject()
2030                rbId = Indx[Obj.GetId()]
2031                if rbId in data['Residue']:
2032                    newNames = []
2033                    newTypes = []
2034                    newXYZ = []
2035                    for i,atype in enumerate(rbData['rbTypes']):
2036                        if atype != 'H':
2037                            newNames.append(rbData['atNames'][i])
2038                            newTypes.append(rbData['rbTypes'][i])
2039                            newXYZ.append(rbData['rbXYZ'][i])
2040                    rbData['atNames'] = newNames
2041                    rbData['rbTypes'] = newTypes
2042                    rbData['rbXYZ'] = newXYZ
2043                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
2044                wx.CallAfter(UpdateResidueRB,rbId)
2045                   
2046            def OnPlotRB(event):
2047                Obj = event.GetEventObject()
2048                Obj.SetValue(False)
2049                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
2050           
2051            nameSizer = wx.BoxSizer(wx.HORIZONTAL)
2052            nameSizer.Add(wx.StaticText(ResidueRBDisplay,-1,'Residue name: '),
2053                0,wx.ALIGN_CENTER_VERTICAL)
2054            RBname = wx.TextCtrl(ResidueRBDisplay,-1,rbData['RBname'])
2055            Indx[RBname.GetId()] = rbId
2056            RBname.Bind(wx.EVT_TEXT_ENTER,OnRBName)
2057            RBname.Bind(wx.EVT_KILL_FOCUS,OnRBName)
2058            nameSizer.Add(RBname,0,wx.ALIGN_CENTER_VERTICAL)
2059            nameSizer.Add((5,0),)
2060            plotRB = wx.CheckBox(ResidueRBDisplay,-1,'Plot?')
2061            Indx[plotRB.GetId()] = rbId
2062            plotRB.Bind(wx.EVT_CHECKBOX,OnPlotRB)
2063            nameSizer.Add(plotRB,0,wx.ALIGN_CENTER_VERTICAL)
2064            nameSizer.Add((5,0),)
2065            if not rbData['useCount']:
2066                delRB = wx.CheckBox(ResidueRBDisplay,-1,'Delete?')
2067                Indx[delRB.GetId()] = rbId
2068                delRB.Bind(wx.EVT_CHECKBOX,OnDelRB)
2069                nameSizer.Add(delRB,0,wx.ALIGN_CENTER_VERTICAL)
2070                if 'H'  in rbData['rbTypes']:
2071                    stripH = wx.CheckBox(ResidueRBDisplay,-1,'Strip H-atoms?')
2072                    Indx[stripH.GetId()] = rbId
2073                    stripH.Bind(wx.EVT_CHECKBOX,OnStripH)
2074                    nameSizer.Add(stripH,0,wx.ALIGN_CENTER_VERTICAL)
2075            return nameSizer
2076           
2077        def rbResidues(rbId,rbData):
2078           
2079            def TypeSelect(event):
2080                AtInfo = data['Residue']['AtInfo']
2081                r,c = event.GetRow(),event.GetCol()
2082                if vecGrid.GetColLabelValue(c) == 'Type':
2083                    PE = G2elemGUI.PickElement(G2frame,oneOnly=True)
2084                    if PE.ShowModal() == wx.ID_OK:
2085                        if PE.Elem != 'None':
2086                            El = PE.Elem.strip().lower().capitalize()
2087                            if El not in AtInfo:
2088                                Info = G2elem.GetAtomInfo(El)
2089                                AtInfo[El] = [Info['Drad']['Color']]
2090                            rbData['rbTypes'][r] = El
2091                            vecGrid.SetCellValue(r,c,El)
2092                    PE.Destroy()
2093
2094            def ChangeCell(event):
2095                r,c =  event.GetRow(),event.GetCol()
2096                if r >= 0 and (0 <= c < 3):
2097                    try:
2098                        val = float(vecGrid.GetCellValue(r,c))
2099                        rbData['rbXYZ'][r][c] = val
2100                    except ValueError:
2101                        pass
2102                       
2103            def RowSelect(event):
2104                r,c =  event.GetRow(),event.GetCol()
2105                if c < 0:                   #only row clicks
2106                    for vecgrid in resList:
2107                        vecgrid.ClearSelection()
2108                    vecGrid.SelectRow(r,True)
2109
2110            def OnRefSel(event):
2111               
2112                Obj = event.GetEventObject()
2113                iref,res,jref = Indx[Obj.GetId()]
2114                sel = Obj.GetValue()
2115                ind = atNames.index(sel)
2116                if rbData['rbTypes'][ind] == 'H':
2117                    G2G.G2MessageBox(G2frame,'You should not select an H-atom for rigid body orientation')
2118                rbData['rbRef'][iref] = ind
2119                FillRefChoice(rbId,rbData)
2120                for i,ref in enumerate(RefObjs[jref]):
2121                    ref.SetItems([atNames[j] for j in refChoice[rbId][i]])
2122                    ref.SetValue(atNames[rbData['rbRef'][i]])                   
2123                rbXYZ = rbData['rbXYZ']
2124                if not iref:     #origin change
2125                    rbXYZ -= rbXYZ[ind]
2126                #TODO - transform all atom XYZ by axis choices
2127                Xxyz = rbXYZ[rbData['rbRef'][1]]
2128                X = Xxyz/np.sqrt(np.sum(Xxyz**2))
2129                Yxyz = rbXYZ[rbData['rbRef'][2]]
2130                Y = Yxyz/np.sqrt(np.sum(Yxyz**2))
2131                Mat = G2mth.getRBTransMat(X,Y)
2132                rbXYZ = np.inner(Mat,rbXYZ).T
2133                rbData['rbXYZ'] = rbXYZ
2134                res.ClearSelection()
2135                resTable = res.GetTable()
2136                for r in range(res.GetNumberRows()):
2137                    row = resTable.GetRowValues(r)
2138                    row[2:4] = rbXYZ[r]
2139                    resTable.SetRowValues(r,row)
2140                res.ForceRefresh()
2141                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
2142               
2143            Types = 2*[wg.GRID_VALUE_STRING,]+3*[wg.GRID_VALUE_FLOAT+':10,5',]
2144            colLabels = ['Name','Type','Cart x','Cart y','Cart z']
2145            table = []
2146            rowLabels = []
2147            for ivec,xyz in enumerate(rbData['rbXYZ']):
2148                table.append([rbData['atNames'][ivec],]+[rbData['rbTypes'][ivec],]+list(xyz))
2149                rowLabels.append(str(ivec))
2150            vecTable = G2G.Table(table,rowLabels=rowLabels,colLabels=colLabels,types=Types)
2151            vecGrid = G2G.GSGrid(ResidueRBDisplay)
2152            Indx[vecGrid.GetId()] = rbId
2153            resList.append(vecGrid)
2154            vecGrid.SetTable(vecTable, True)
2155            if 'phoenix' in wx.version():
2156                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGED, ChangeCell)
2157            else:
2158                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGE, ChangeCell)
2159            vecGrid.Bind(wg.EVT_GRID_CELL_LEFT_DCLICK, TypeSelect)
2160            vecGrid.Bind(wg.EVT_GRID_LABEL_LEFT_CLICK, RowSelect)
2161            attr = wx.grid.GridCellAttr()
2162            attr.SetEditor(G2G.GridFractionEditor(vecGrid))
2163            for c in range(3):
2164                vecGrid.SetColAttr(c, attr)
2165            for row in range(vecTable.GetNumberRows()):
2166                for col in range(5):
2167                    vecGrid.SetCellStyle(row,col,VERY_LIGHT_GREY,True)
2168            vecGrid.AutoSizeColumns(False)
2169            vecSizer = wx.BoxSizer()
2170            vecSizer.Add(vecGrid)
2171           
2172            refAtmSizer = wx.BoxSizer(wx.HORIZONTAL)
2173            atNames = rbData['atNames']
2174            rbRef = rbData['rbRef']
2175            if rbData['rbRef'][3] or rbData['useCount']:
2176                refAtmSizer.Add(wx.StaticText(ResidueRBDisplay,-1,
2177                    'Orientation reference non-H atoms A-B-C: %s, %s, %s'%(atNames[rbRef[0]], \
2178                     atNames[rbRef[1]],atNames[rbRef[2]])),0)
2179            else:
2180                refAtmSizer.Add(wx.StaticText(ResidueRBDisplay,-1,
2181                    'Orientation reference non-H atoms A-B-C: '),0,wx.ALIGN_CENTER_VERTICAL)
2182                refObj = [0,0,0]
2183                for i in range(3):
2184                    choices = [atNames[j] for j in refChoice[rbId][i]]
2185                    refSel = wx.ComboBox(ResidueRBDisplay,-1,value='',
2186                        choices=choices,style=wx.CB_READONLY|wx.CB_DROPDOWN)
2187                    refSel.SetValue(atNames[rbRef[i]])
2188                    refSel.Bind(wx.EVT_COMBOBOX, OnRefSel)
2189                    Indx[refSel.GetId()] = [i,vecGrid,len(RefObjs)]
2190                    refObj[i] = refSel
2191                    refAtmSizer.Add(refSel,0,wx.ALIGN_CENTER_VERTICAL)
2192                RefObjs.append(refObj)
2193           
2194            mainSizer = wx.BoxSizer(wx.VERTICAL)
2195            mainSizer.Add(refAtmSizer)
2196            mainSizer.Add(vecSizer)
2197            return mainSizer
2198           
2199        def SeqSizer(angSlide,rbId,iSeq,Seq,atNames):
2200           
2201            def ChangeAngle(event):
2202                event.Skip()
2203                Obj = event.GetEventObject()
2204                rbId,Seq = Indx[Obj.GetId()][:2]
2205                val = Seq[2]
2206                try:
2207                    val = float(Obj.GetValue())
2208                    Seq[2] = val
2209                except ValueError:
2210                    pass
2211                Obj.SetValue('%8.2f'%(val))
2212                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,data['Residue'][rbId],plotDefaults)
2213               
2214            def OnRadBtn(event):
2215                Obj = event.GetEventObject()
2216                Seq,iSeq,angId = Indx[Obj.GetId()]
2217                data['Residue'][rbId]['SelSeq'] = [iSeq,angId]
2218                angSlide.SetValue(int(100*Seq[2]))
2219               
2220            def OnDelBtn(event):
2221                Obj = event.GetEventObject()
2222                rbId,Seq = Indx[Obj.GetId()]
2223                data['Residue'][rbId]['rbSeq'].remove(Seq)       
2224                wx.CallAfter(UpdateResidueRB,rbId)
2225           
2226            seqSizer = wx.FlexGridSizer(0,5,2,2)
2227            seqSizer.AddGrowableCol(3,0)
2228            iBeg,iFin,angle,iMove = Seq
2229            ang = wx.TextCtrl(ResidueRBDisplay,-1,'%8.2f'%(angle),size=(50,20))
2230            if not iSeq:
2231                radBt = wx.RadioButton(ResidueRBDisplay,-1,'',style=wx.RB_GROUP)
2232                data['Residue'][rbId]['SelSeq'] = [iSeq,ang.GetId()]
2233            else:
2234                radBt = wx.RadioButton(ResidueRBDisplay,-1,'')
2235            radBt.Bind(wx.EVT_RADIOBUTTON,OnRadBtn)                   
2236            seqSizer.Add(radBt)
2237            delBt = wx.RadioButton(ResidueRBDisplay,-1,'')
2238            delBt.Bind(wx.EVT_RADIOBUTTON,OnDelBtn)
2239            seqSizer.Add(delBt)
2240            bond = wx.TextCtrl(ResidueRBDisplay,-1,'%s %s'%(atNames[iBeg],atNames[iFin]),size=(50,20))
2241            seqSizer.Add(bond,0,wx.ALIGN_CENTER_VERTICAL)
2242            Indx[radBt.GetId()] = [Seq,iSeq,ang.GetId()]
2243            Indx[delBt.GetId()] = [rbId,Seq]
2244            Indx[ang.GetId()] = [rbId,Seq,ang]
2245            ang.Bind(wx.EVT_TEXT_ENTER,ChangeAngle)
2246            ang.Bind(wx.EVT_KILL_FOCUS,ChangeAngle)
2247            seqSizer.Add(ang,0,wx.ALIGN_CENTER_VERTICAL)
2248            atms = ''
2249            for i in iMove:   
2250                atms += ' %s,'%(atNames[i])
2251            moves = wx.TextCtrl(ResidueRBDisplay,-1,atms[:-1],size=(200,20))
2252            seqSizer.Add(moves,1,wx.ALIGN_CENTER_VERTICAL|wx.EXPAND|wx.RIGHT)
2253            return seqSizer
2254           
2255        def SlideSizer():
2256           
2257            def OnSlider(event):
2258                Obj = event.GetEventObject()
2259                rbData = Indx[Obj.GetId()]
2260                iSeq,angId = rbData['SelSeq']
2261                val = float(Obj.GetValue())/100.
2262                rbData['rbSeq'][iSeq][2] = val
2263                Indx[angId][2].SetValue('%8.2f'%(val))
2264                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
2265           
2266            slideSizer = wx.BoxSizer(wx.HORIZONTAL)
2267            slideSizer.Add(wx.StaticText(ResidueRBDisplay,-1,'Selected torsion angle:'),0)
2268            iSeq,angId = rbData['SelSeq']
2269            angSlide = wx.Slider(ResidueRBDisplay,-1,
2270                int(100*rbData['rbSeq'][iSeq][2]),0,36000,size=(200,20),
2271                style=wx.SL_HORIZONTAL)
2272            angSlide.Bind(wx.EVT_SLIDER, OnSlider)
2273            Indx[angSlide.GetId()] = rbData
2274            slideSizer.Add(angSlide,0)           
2275            return slideSizer,angSlide
2276           
2277        def FillRefChoice(rbId,rbData):
2278            choiceIds = [i for i in range(len(rbData['atNames']))]
2279            for seq in rbData['rbSeq']:
2280                for i in seq[3]:
2281                    try:
2282                        choiceIds.remove(i)
2283                    except ValueError:
2284                        pass
2285            rbRef = rbData['rbRef']
2286            for i in range(3):
2287                try:
2288                    choiceIds.remove(rbRef[i])
2289                except ValueError:
2290                    pass
2291            refChoice[rbId] = [choiceIds[:],choiceIds[:],choiceIds[:]]
2292            for i in range(3):
2293                refChoice[rbId][i].append(rbRef[i])
2294                refChoice[rbId][i].sort()
2295               
2296        def OnSelect(event):
2297            rbname = rbchoice[select.GetSelection()]
2298            rbId = RBnames[rbname]
2299            wx.CallLater(100,UpdateResidueRB,rbId)
2300           
2301        GS = ResidueRBDisplay.GetSizer()
2302        if GS: 
2303            try:        #get around a c++ error in wx 4.0; doing is again seems to be OK
2304                GS.Clear(True)
2305            except:
2306                GS.Clear(True)
2307       
2308        RBnames = {}
2309        for rbid in data['RBIds']['Residue']:
2310            RBnames.update({data['Residue'][rbid]['RBname']:rbid,})
2311        if not RBnames:
2312            return
2313        rbchoice = list(RBnames.keys())
2314        ResidueRBSizer = wx.BoxSizer(wx.VERTICAL)
2315        if len(RBnames) > 1:
2316            selSizer = wx.BoxSizer(wx.HORIZONTAL)
2317            selSizer.Add(wx.StaticText(ResidueRBDisplay,label=' Select residue to view:'),0)
2318            rbchoice.sort()
2319            select = wx.ComboBox(ResidueRBDisplay,choices=rbchoice)
2320            select.Bind(wx.EVT_COMBOBOX,OnSelect)
2321            selSizer.Add(select,0)
2322            ResidueRBSizer.Add(selSizer,0)
2323        if not rbId:
2324            rbId = RBnames[rbchoice[0]]
2325        rbData = data['Residue'][rbId]
2326        FillRefChoice(rbId,rbData)
2327        ResidueRBSizer.Add(rbNameSizer(rbId,rbData),0)
2328        ResidueRBSizer.Add(rbResidues(rbId,rbData),0)
2329        ResidueRBSizer.Add((5,5),0)
2330        if rbData['rbSeq']:
2331            slideSizer,angSlide = SlideSizer()
2332        if len(rbData['rbSeq']):
2333            ResidueRBSizer.Add(wx.StaticText(ResidueRBDisplay,-1,
2334                'Sel  Del  Bond             Angle      Riding atoms'),
2335                0,wx.ALIGN_CENTER_VERTICAL)                       
2336        for iSeq,Seq in enumerate(rbData['rbSeq']):
2337            ResidueRBSizer.Add(SeqSizer(angSlide,rbId,iSeq,Seq,rbData['atNames']))
2338        if rbData['rbSeq']:
2339            ResidueRBSizer.Add(slideSizer,)
2340
2341        ResidueRBSizer.Add((5,25),)
2342        ResidueRBSizer.Layout()   
2343        ResidueRBDisplay.SetSizer(ResidueRBSizer,True)
2344        ResidueRBDisplay.SetAutoLayout(True)
2345        Size = ResidueRBSizer.GetMinSize()
2346        ResidueRBDisplay.SetSize(Size)
2347        ResidueRBDisplay.Show()
2348       
2349    def SetStatusLine(text):
2350        G2frame.GetStatusBar().SetStatusText(text,1)                                     
2351
2352    G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.RigidBodyMenu)
2353    SetStatusLine('')
2354    UpdateVectorRB()
2355    G2frame.rbBook.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, OnPageChanged)
2356    wx.CallAfter(OnPageChanged,None)
Note: See TracBrowser for help on using the repository browser.