source: trunk/GSASIIconstrGUI.py @ 3944

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

fix for vector rigid body display

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 103.9 KB
Line 
1# -*- coding: utf-8 -*-
2#GSASIIconstrGUI - constraint GUI routines
3########### SVN repository information ###################
4# $Date: 2019-05-02 18:00:36 +0000 (Thu, 02 May 2019) $
5# $Author: toby $
6# $Revision: 3944 $
7# $URL: trunk/GSASIIconstrGUI.py $
8# $Id: GSASIIconstrGUI.py 3944 2019-05-02 18:00:36Z toby $
9########### SVN repository information ###################
10'''
11*GSASIIconstrGUI: Constraint GUI routines*
12------------------------------------------
13
14Used to define constraints and rigid bodies.
15
16'''
17from __future__ import division, print_function
18import 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: 3944 $")
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    invTrans = nl.inv(Trans)
1339    opId = oldPhase['pId']
1340    npId = newPhase['pId']
1341    cx,ct,cs,cia = newPhase['General']['AtomPtrs']
1342    nAtoms = newPhase['Atoms']
1343    oSGData = oldPhase['General']['SGData']
1344    nSGData = newPhase['General']['SGData']
1345    #oAcof = G2lat.cell2A(oldPhase['General']['Cell'][1:7])
1346    #nAcof = G2lat.cell2A(newPhase['General']['Cell'][1:7])
1347    item = G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Constraints')
1348    if not item:
1349        print('Error: no constraints in Data Tree')
1350        return
1351    constraints = G2frame.GPXtree.GetItemPyData(item)
1352    xnames = ['dAx','dAy','dAz']
1353    # constraints on matching atom params between phases
1354    for ia,code in enumerate(atCodes):
1355        atom = nAtoms[ia]
1356        if not ia and atom[cia] == 'A':
1357            wx.MessageDialog(G2frame,
1358                'Anisotropic thermal motion constraints are not developed at the present time',
1359                'Anisotropic thermal constraint?',style=wx.ICON_INFORMATION).ShowModal()
1360        siteSym = G2spc.SytSym(atom[cx:cx+3],nSGData)[0]
1361        CSX = G2spc.GetCSxinel(siteSym)
1362#        CSU = G2spc.GetCSuinel(siteSym)
1363        item = code.split('+')[0]
1364        iat,opr = item.split(':')
1365        Nop = abs(int(opr))%100-1
1366        if '-' in opr:
1367            Nop *= -1
1368        Opr = oldPhase['General']['SGData']['SGOps'][abs(Nop)][0]
1369        if Nop < 0:         #inversion
1370            Opr *= -1
1371        XOpr = np.inner(Opr,Trans)
1372        for i,ix in enumerate(list(CSX[0])):
1373            if not ix:
1374                continue
1375            name = xnames[i]
1376            IndpCon = [1.0,G2obj.G2VarObj('%d::%s:%d'%(npId,name,ia))]
1377            DepCons = []
1378            for iop,opval in enumerate(XOpr[i]):
1379                if opval:
1380                    DepCons.append([opval,G2obj.G2VarObj('%d::%s:%s'%(opId,xnames[iop],iat))])
1381            if len(DepCons) == 1:
1382                constraints['Phase'].append([DepCons[0],IndpCon,None,None,'e'])
1383            elif len(DepCons) > 1:
1384                for Dep in DepCons:
1385                    Dep[0] *= -1
1386                constraints['Phase'].append([IndpCon]+DepCons+[0.0,None,'c'])
1387        for name in ['Afrac','AUiso']:
1388            IndpCon = [1.0,G2obj.G2VarObj('%d::%s:%d'%(npId,name,ia))]
1389            DepCons = [1.0,G2obj.G2VarObj('%d::%s:%s'%(opId,name,iat))]
1390            constraints['Phase'].append([DepCons,IndpCon,None,None,'e'])
1391           
1392        # unfinished Anisotropic constraint generation
1393#        Uids = [[0,0,'AU11'],[1,1,'AU22'],[2,2,'AU33'],[0,1,'AU12'],[0,2,'AU13'],[1,2,'AU23']]
1394#        DepConsDict = dict(zip(Us,[[],[],[],[],[],[]]))
1395#        for iu,Uid in enumerate(Uids):
1396#            UMT = np.zeros((3,3))
1397#            UMT[Uid[0],Uid[1]] = 1
1398#            nUMT = G2lat.prodMGMT(UMT,invTrans)
1399#            nUT = G2lat.UijtoU6(nUMT)
1400#            for iu,nU in enumerate(nUT):
1401#                if abs(nU) > 1.e-8:
1402#                    parm = '%d::%s;%s'%(opId,Us[iu],iat)
1403#                    DepConsDict[Uid[2]].append([abs(nU%1.),G2obj.G2VarObj(parm)])
1404#        nUcof = atom[iu:iu+6]
1405#        conStrings = []
1406#        for iU,Usi in enumerate(Us):
1407#            parm = '%d::%s;%d'%(npId,Usi,ia)
1408#            parmDict[parm] = nUcof[iU]
1409#            varyList.append(parm)
1410#            IndpCon = [1.0,G2obj.G2VarObj(parm)]
1411#            conStr = str([IndpCon,DepConsDict[Usi]])
1412#            if conStr in conStrings:
1413#                continue
1414#            conStrings.append(conStr)
1415#            if len(DepConsDict[Usi]) == 1:
1416#                if DepConsDict[Usi][0]:
1417#                    constraints['Phase'].append([IndpCon,DepConsDict[Usi][0],None,None,'e'])
1418#            elif len(DepConsDict[Usi]) > 1:       
1419#                for Dep in DepConsDict[Usi]:
1420#                    Dep[0] *= -1
1421#                constraints['Phase'].append([IndpCon]+DepConsDict[Usi]+[0.0,None,'c'])
1422           
1423        #how do I do Uij's for most Trans?
1424
1425    # constraints on lattice parameters between phases
1426#    T = nl.inv(Trans).T
1427    T = Trans.T
1428    conMat = [
1429        [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]],
1430        [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]],
1431        [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]],
1432        [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]],
1433        [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]],
1434        [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]]
1435        ]
1436    # Gnew = conMat * A:
1437#         T00**2*a0  T01**2*a1 T02**2*a2 T00*T01*a3    T00*T02*a4    T01*T02*a5
1438#         T10**2*a0  T11**2*a1 T12**2*a2 T10*T11*a3    T10*T12*a4    T11*T12*a5
1439#         T20**2*a0  T21**2*a1 T22**2*a2 T20*T21*a3    T20*T22*a4    T21*T22*a5
1440#         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
1441#         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
1442#         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
1443    # Generated as symbolic code using:
1444    # import sym
1445    # A0, A1, A2, A3, A4, A5 = sym.symbols('A0, A1, A2, A3, A4, A5')
1446    # G = sym.Matrix([ [A0,    A3/2.,  A4/2.], [A3/2.,  A1,    A5/2.], [A4/2., A5/2.,    A2]])
1447    # transformation matrix
1448    # T00, T10, T20, T01, T11, T21, T02, T12, T22 = sym.symbols('T00, T10, T20, T01, T11, T21, T02, T12, T22')
1449    # Tr = sym.Matrix([ [T00, T10, T20], [T01, T11, T21], [T02, T12, T22],])
1450    # Gnew = (Tr.T*G)*Tr
1451   
1452    #print('old A',G2lat.cell2A(oldPhase['General']['Cell'][1:7]))
1453    #print('new A',G2lat.cell2A(newPhase['General']['Cell'][1:7]))
1454    for iAnew,Asi in enumerate(['A0','A1','A2','A3','A4','A5']): # loop through A[i] for new cell
1455        Nparm = str(npId) + '::' + Asi
1456        if Nparm != SetUniqAj(npId,iAnew,nSGData):
1457            continue # skip: Ai constrained from Aj or must be zero
1458        multDict = {}
1459        for iAorg in range(6):
1460            cA = conMat[iAnew][iAorg] # coeff for A[i] in constraint matrix
1461            if abs(cA) < 1.e-8: continue
1462            parm = SetUniqAj(opId,iAorg,oSGData) # translate to unique A[i] in original cell
1463            if not parm: continue # must be zero
1464            # sum coeff
1465            if parm in multDict:
1466                multDict[parm] += cA
1467            else:
1468                multDict[parm] = cA
1469        # any non-zero multipliers?
1470        maxMult = 0
1471        for i in multDict:
1472            maxMult = max(maxMult,abs(multDict[i]))
1473        if maxMult <= 0:  # Nparm computes as zero; Fix this parameter
1474            constraints['Phase'] += [[
1475                [0.0,G2obj.G2VarObj(Nparm)],
1476                None,None,'h']]
1477        elif len(multDict) == 1:        # create equivalence
1478            key = list(multDict.keys())[0]
1479            constraints['Phase'] += [[
1480                [1.0,G2obj.G2VarObj(key)],
1481                [multDict[key],G2obj.G2VarObj(Nparm)],
1482                None,None,'e']]
1483        else:                           # create constraint
1484            constr = [[-1.0,G2obj.G2VarObj(Nparm)]]
1485            for key in multDict:
1486                constr += [[multDict[key],G2obj.G2VarObj(key)]]
1487            constr += [0.0,None,'c']
1488            constraints['Phase'] += [constr]
1489   
1490    # constraints on HAP Scale, etc.
1491    for hId,hist in enumerate(UseList):    #HAP - seems OK
1492        ohapkey = '%d:%d:'%(opId,hId)
1493        nhapkey = '%d:%d:'%(npId,hId)
1494        IndpCon = [1.0,G2obj.G2VarObj(ohapkey+'Scale')]
1495        DepCons = [detTrans,G2obj.G2VarObj(nhapkey+'Scale')]
1496        constraints['HAP'].append([DepCons,IndpCon,None,None,'e'])
1497        for name in ['Size;i','Mustrain;i']:
1498            IndpCon = [1.0,G2obj.G2VarObj(ohapkey+name)]
1499            DepCons = [1.0,G2obj.G2VarObj(nhapkey+name)]
1500            constraints['HAP'].append([IndpCon,DepCons,None,None,'e'])
1501       
1502################################################################################
1503#### Rigid bodies
1504################################################################################
1505
1506def UpdateRigidBodies(G2frame,data):
1507    '''Called when Rigid bodies tree item is selected.
1508    Displays the rigid bodies in the data window
1509    '''
1510    if not data.get('RBIds') or not data:
1511        data.update({'Vector':{'AtInfo':{}},'Residue':{'AtInfo':{}},
1512            'RBIds':{'Vector':[],'Residue':[]}})       #empty/bad dict - fill it
1513           
1514    global resList,rbId
1515    Indx = {}
1516    resList = []
1517    plotDefaults = {'oldxy':[0.,0.],'Quaternion':[0.,0.,0.,1.],'cameraPos':30.,'viewDir':[0,0,1],}
1518
1519    G2frame.rbBook = G2G.GSNoteBook(parent=G2frame.dataWindow)
1520    G2frame.dataWindow.GetSizer().Add(G2frame.rbBook,1,wx.ALL|wx.EXPAND)
1521    VectorRB = wx.ScrolledWindow(G2frame.rbBook)
1522    VectorRBDisplay = wx.Panel(VectorRB)
1523    G2frame.rbBook.AddPage(VectorRB,'Vector rigid bodies')
1524    ResidueRB = wx.ScrolledWindow(G2frame.rbBook)
1525    ResidueRBDisplay = wx.Panel(ResidueRB)
1526    G2frame.rbBook.AddPage(ResidueRB,'Residue rigid bodies')
1527   
1528    def OnPageChanged(event):
1529        global resList
1530        resList = []
1531        if event:       #page change event!
1532            page = event.GetSelection()
1533        else:
1534            page = G2frame.rbBook.GetSelection()
1535        G2frame.rbBook.ChangeSelection(page)
1536        text = G2frame.rbBook.GetPageText(page)
1537        if text == 'Vector rigid bodies':
1538            G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.VectorBodyMenu)
1539            G2frame.Bind(wx.EVT_MENU, AddVectorRB, id=G2G.wxID_VECTORBODYADD)
1540            G2frame.Page = [page,'vrb']
1541            UpdateVectorRB()
1542        elif text == 'Residue rigid bodies':
1543            G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.RigidBodyMenu)
1544            G2frame.Bind(wx.EVT_MENU, AddResidueRB, id=G2G.wxID_RIGIDBODYADD)
1545            G2frame.Bind(wx.EVT_MENU, OnImportRigidBody, id=G2G.wxID_RIGIDBODYIMPORT)
1546            G2frame.Bind(wx.EVT_MENU, OnDefineTorsSeq, id=G2G.wxID_RESIDUETORSSEQ) #enable only if residue RBs exist?
1547            G2frame.Page = [page,'rrb']
1548            UpdateResidueRB()
1549           
1550    def getMacroFile(macName):
1551        defDir = os.path.join(os.path.split(__file__)[0],'GSASIImacros')
1552        dlg = wx.FileDialog(G2frame,message='Choose '+macName+' rigid body macro file',
1553            defaultDir=defDir,defaultFile="",wildcard="GSAS-II macro file (*.mac)|*.mac",
1554            style=wx.FD_OPEN | wx.FD_CHANGE_DIR)
1555        try:
1556            if dlg.ShowModal() == wx.ID_OK:
1557                macfile = dlg.GetPath()
1558                macro = open(macfile,'Ur')
1559                head = macro.readline()
1560                if macName not in head:
1561                    print (head)
1562                    print ('**** ERROR - wrong restraint macro file selected, try again ****')
1563                    macro = []
1564            else: # cancel was pressed
1565                macro = []
1566        finally:
1567            dlg.Destroy()
1568        return macro        #advanced past 1st line
1569       
1570    def getTextFile():
1571        dlg = wx.FileDialog(G2frame,'Choose rigid body text file', '.', '',
1572            "GSAS-II text file (*.txt)|*.txt|XYZ file (*.xyz)|*.xyz|"
1573            "Sybyl mol2 file (*.mol2)|*.mol2|PDB file (*.pdb;*.ent)|*.pdb;*.ent",
1574            wx.FD_OPEN | wx.FD_CHANGE_DIR)
1575        try:
1576            if dlg.ShowModal() == wx.ID_OK:
1577                txtfile = dlg.GetPath()
1578                ext = os.path.splitext(txtfile)[1]
1579                text = open(txtfile,'Ur')
1580            else: # cancel was pressed
1581                ext = ''
1582                text = []
1583        finally:
1584            dlg.Destroy()
1585        if 'ent' in ext:
1586            ext = '.pdb'
1587        return text,ext.lower()
1588       
1589    def OnImportRigidBody(event):
1590        page = G2frame.rbBook.GetSelection()
1591        if 'Vector' in G2frame.rbBook.GetPageText(page):
1592            pass
1593        elif 'Residue' in G2frame.rbBook.GetPageText(page):
1594            ImportResidueRB()
1595           
1596    def AddVectorRB(event):
1597        AtInfo = data['Vector']['AtInfo']
1598        dlg = G2G.MultiIntegerDialog(G2frame,'New Rigid Body',['No. atoms','No. translations'],[1,1])
1599        if dlg.ShowModal() == wx.ID_OK:
1600            nAtoms,nTrans = dlg.GetValues()
1601            rbId = ran.randint(0,sys.maxsize)
1602            vecMag = [1.0 for i in range(nTrans)]
1603            vecRef = [False for i in range(nTrans)]
1604            vecVal = [np.zeros((nAtoms,3)) for j in range(nTrans)]
1605            rbTypes = ['C' for i in range(nAtoms)]
1606            Info = G2elem.GetAtomInfo('C')
1607            AtInfo['C'] = [Info['Drad'],Info['Color']]
1608            data['Vector'][rbId] = {'RBname':'UNKRB','VectMag':vecMag,'rbXYZ':np.zeros((nAtoms,3)),
1609                'rbRef':[0,1,2,False],'VectRef':vecRef,'rbTypes':rbTypes,'rbVect':vecVal,'useCount':0}
1610            data['RBIds']['Vector'].append(rbId)
1611        dlg.Destroy()
1612        UpdateVectorRB()
1613       
1614    def AddResidueRB(event):
1615        global rbId
1616        AtInfo = data['Residue']['AtInfo']
1617        macro = getMacroFile('rigid body')
1618        if not macro:
1619            return
1620        macStr = macro.readline()
1621        while macStr:
1622            items = macStr.split()
1623            if 'I' == items[0]:
1624                rbId = ran.randint(0,sys.maxsize)
1625                rbName = items[1]
1626                rbTypes = []
1627                rbXYZ = []
1628                rbSeq = []
1629                atNames = []
1630                nAtms,nSeq,nOrig,mRef,nRef = [int(items[i]) for i in [2,3,4,5,6]]
1631                for iAtm in range(nAtms):
1632                    macStr = macro.readline().split()
1633                    atName = macStr[0]
1634                    atType = macStr[1]
1635                    atNames.append(atName)
1636                    rbXYZ.append([float(macStr[i]) for i in [2,3,4]])
1637                    rbTypes.append(atType)
1638                    if atType not in AtInfo:
1639                        Info = G2elem.GetAtomInfo(atType)
1640                        AtInfo[atType] = [Info['Drad'],Info['Color']]
1641                rbXYZ = np.array(rbXYZ)-np.array(rbXYZ[nOrig-1])
1642                for iSeq in range(nSeq):
1643                    macStr = macro.readline().split()
1644                    mSeq = int(macStr[0])
1645                    for jSeq in range(mSeq):
1646                        macStr = macro.readline().split()
1647                        iBeg = int(macStr[0])-1
1648                        iFin = int(macStr[1])-1
1649                        angle = 0.0
1650                        nMove = int(macStr[2])
1651                        iMove = [int(macStr[i])-1 for i in range(3,nMove+3)]
1652                        rbSeq.append([iBeg,iFin,angle,iMove])
1653                data['Residue'][rbId] = {'RBname':rbName,'rbXYZ':rbXYZ,'rbTypes':rbTypes,
1654                    'atNames':atNames,'rbRef':[nOrig-1,mRef-1,nRef-1,True],'rbSeq':rbSeq,
1655                    'SelSeq':[0,0],'useCount':0}
1656                data['RBIds']['Residue'].append(rbId)
1657                print ('Rigid body '+rbName+' added')
1658            macStr = macro.readline()
1659        macro.close()
1660        UpdateResidueRB()
1661       
1662    def ImportResidueRB():
1663        global rbId
1664        AtInfo = data['Residue']['AtInfo']
1665        text,ext = getTextFile()
1666        if not text:
1667            return
1668        rbId = ran.randint(0,sys.maxsize)
1669        rbTypes = []
1670        rbXYZ = []
1671        atNames = []
1672        txtStr = text.readline()
1673        if 'xyz' in ext:
1674            txtStr = text.readline()
1675            txtStr = text.readline()
1676        elif 'mol2' in ext:
1677            while 'ATOM' not in txtStr:
1678                txtStr = text.readline()
1679            txtStr = text.readline()
1680        elif 'pdb' in ext:
1681            while 'ATOM' not in txtStr[:6] and 'HETATM' not in txtStr[:6]:
1682                txtStr = text.readline()
1683                #print txtStr
1684        items = txtStr.split()
1685        while len(items):
1686            if 'txt' in ext:
1687                atName = items[0]
1688                atType = items[1]
1689                rbXYZ.append([float(items[i]) for i in [2,3,4]])
1690            elif 'xyz' in ext:
1691                atType = items[0]
1692                rbXYZ.append([float(items[i]) for i in [1,2,3]])
1693                atName = atType+str(len(rbXYZ))
1694            elif 'mol2' in ext:
1695                atType = items[1]
1696                atName = items[1]+items[0]
1697                rbXYZ.append([float(items[i]) for i in [2,3,4]])
1698            elif 'pdb' in ext:
1699                atType = items[-1]
1700                atName = items[2]
1701                xyz = txtStr[30:55].split()                   
1702                rbXYZ.append([float(x) for x in xyz])
1703            atNames.append(atName)
1704            rbTypes.append(atType)
1705            if atType not in AtInfo:
1706                Info = G2elem.GetAtomInfo(atType)
1707                AtInfo[atType] = [Info['Drad'],Info['Color']]
1708            txtStr = text.readline()
1709            if 'mol2' in ext and 'BOND' in txtStr:
1710                break
1711            if 'pdb' in ext and ('ATOM' not in txtStr[:6] and 'HETATM' not in txtStr[:6]):
1712                break
1713            items = txtStr.split()
1714        if len(atNames) < 3:
1715            G2G.G2MessageBox(G2frame,'Not enough atoms in rigid body; must be 3 or more')
1716        else:
1717            rbXYZ = np.array(rbXYZ)-np.array(rbXYZ[0])
1718            Xxyz = rbXYZ[1]
1719            X = Xxyz/np.sqrt(np.sum(Xxyz**2))
1720            Yxyz = rbXYZ[2]
1721            Y = Yxyz/np.sqrt(np.sum(Yxyz**2))
1722            Mat = G2mth.getRBTransMat(X,Y)
1723            rbXYZ = np.inner(Mat,rbXYZ).T
1724            data['Residue'][rbId] = {'RBname':'UNKRB','rbXYZ':rbXYZ,'rbTypes':rbTypes,
1725                'atNames':atNames,'rbRef':[0,1,2,False],'rbSeq':[],'SelSeq':[0,0],'useCount':0}
1726            data['RBIds']['Residue'].append(rbId)
1727            print ('Rigid body UNKRB added')
1728        text.close()
1729        UpdateResidueRB(rbId)
1730       
1731    def FindNeighbors(Orig,XYZ,atTypes,atNames,AtInfo):
1732        Radii = []
1733        for Atype in atTypes:
1734            Radii.append(AtInfo[Atype][0])
1735        Radii = np.array(Radii)
1736        Neigh = []
1737        Dx = XYZ-XYZ[Orig]
1738        dist = np.sqrt(np.sum(Dx**2,axis=1))
1739        sumR = Radii[Orig]+Radii
1740        IndB = ma.nonzero(ma.masked_greater(dist-0.85*sumR,0.))
1741        for j in IndB[0]:
1742            if j != Orig and atTypes[j] != 'H':
1743                Neigh.append(atNames[j])
1744        return Neigh
1745       
1746    def FindAllNeighbors(XYZ,atTypes,atNames,AtInfo):
1747        NeighDict = {}
1748        for iat,xyz in enumerate(atNames):
1749            NeighDict[atNames[iat]] = FindNeighbors(iat,XYZ,atTypes,atNames,AtInfo)
1750        return NeighDict
1751       
1752    def FindRiding(Orig,Pivot,NeighDict):
1753        riding = [Orig,Pivot]
1754        iAdd = 1
1755        new = True
1756        while new:
1757            newAtms = NeighDict[riding[iAdd]]
1758            for At in newAtms:
1759                new = False
1760                if At not in riding:
1761                    riding.append(At)
1762                    new = True
1763            iAdd += 1
1764            if iAdd < len(riding):
1765                new = True
1766        return riding[2:]
1767                       
1768    def OnDefineTorsSeq(event):
1769        global rbId
1770        rbData = data['Residue'][rbId]
1771        if not len(rbData):
1772            return
1773        atNames = rbData['atNames']
1774        AtInfo = data['Residue']['AtInfo']
1775        atTypes = rbData['rbTypes']
1776        XYZ = rbData['rbXYZ']
1777        neighDict = FindAllNeighbors(XYZ,atTypes,atNames,AtInfo)
1778        TargList = []           
1779        dlg = wx.SingleChoiceDialog(G2frame,'Select origin atom for torsion sequence','Origin atom',rbData['atNames'])
1780        if dlg.ShowModal() == wx.ID_OK:
1781            Orig = dlg.GetSelection()
1782            TargList = neighDict[atNames[Orig]]
1783        dlg.Destroy()
1784        if not len(TargList):
1785            return
1786        dlg = wx.SingleChoiceDialog(G2frame,'Select pivot atom for torsion sequence','Pivot atom',TargList)
1787        if dlg.ShowModal() == wx.ID_OK:
1788            Piv = atNames.index(TargList[dlg.GetSelection()])
1789            riding = FindRiding(atNames[Orig],atNames[Piv],neighDict)
1790            Riding = []
1791            for atm in riding:
1792                Riding.append(atNames.index(atm))
1793            rbData['rbSeq'].append([Orig,Piv,0.0,Riding])           
1794        dlg.Destroy()
1795        UpdateResidueRB(rbId)
1796
1797    def UpdateVectorRB(Scroll=0):
1798        AtInfo = data['Vector']['AtInfo']
1799        refChoice = {}
1800        if 'DELETED' in str(G2frame.GetStatusBar()):   #seems to be no other way to do this (wx bug)
1801            if GSASIIpath.GetConfigValue('debug'):
1802                print ('DBG_wx error: Rigid Body/Status not cleanly deleted after Refine')
1803            return
1804        SetStatusLine(' You may use e.g. "c60" or "s60" for a vector entry')
1805        def rbNameSizer(rbId,rbData):
1806
1807            def OnRBName(event):
1808                event.Skip()
1809                Obj = event.GetEventObject()
1810                rbData['RBname'] = Obj.GetValue()
1811               
1812            def OnDelRB(event):
1813                Obj = event.GetEventObject()
1814                rbId = Indx[Obj.GetId()]
1815                if rbId in data['Vector']:
1816                    del data['Vector'][rbId]
1817                    data['RBIds']['Vector'].remove(rbId)
1818                    rbData['useCount'] -= 1
1819                wx.CallAfter(UpdateVectorRB)
1820               
1821            def OnPlotRB(event):
1822                Obj = event.GetEventObject()
1823                Obj.SetValue(False)
1824                G2plt.PlotRigidBody(G2frame,'Vector',AtInfo,rbData,plotDefaults)
1825           
1826            nameSizer = wx.BoxSizer(wx.HORIZONTAL)
1827            nameSizer.Add(wx.StaticText(VectorRBDisplay,-1,'Rigid body name: '),
1828                0,wx.ALIGN_CENTER_VERTICAL)
1829            RBname = wx.TextCtrl(VectorRBDisplay,-1,rbData['RBname'])
1830            Indx[RBname.GetId()] = rbId
1831            RBname.Bind(wx.EVT_TEXT_ENTER,OnRBName)
1832            RBname.Bind(wx.EVT_KILL_FOCUS,OnRBName)
1833            nameSizer.Add(RBname,0,wx.ALIGN_CENTER_VERTICAL)
1834            nameSizer.Add((5,0),)
1835            plotRB = wx.CheckBox(VectorRBDisplay,-1,'Plot?')
1836            Indx[plotRB.GetId()] = rbId
1837            plotRB.Bind(wx.EVT_CHECKBOX,OnPlotRB)
1838            nameSizer.Add(plotRB,0,wx.ALIGN_CENTER_VERTICAL)
1839            nameSizer.Add((5,0),)
1840            if not rbData['useCount']:
1841                delRB = wx.CheckBox(VectorRBDisplay,-1,'Delete?')
1842                Indx[delRB.GetId()] = rbId
1843                delRB.Bind(wx.EVT_CHECKBOX,OnDelRB)
1844                nameSizer.Add(delRB,0,wx.ALIGN_CENTER_VERTICAL)
1845            return nameSizer
1846           
1847        def rbRefAtmSizer(rbId,rbData):
1848           
1849            def OnRefSel(event):
1850                Obj = event.GetEventObject()
1851                iref = Indx[Obj.GetId()]
1852                sel = Obj.GetValue()
1853                rbData['rbRef'][iref] = atNames.index(sel)
1854                FillRefChoice(rbId,rbData)
1855           
1856            refAtmSizer = wx.BoxSizer(wx.HORIZONTAL)
1857            atNames = [name+str(i) for i,name in enumerate(rbData['rbTypes'])]
1858            rbRef = rbData.get('rbRef',[0,1,2,False])
1859            rbData['rbRef'] = rbRef
1860            if rbData['useCount']:
1861                refAtmSizer.Add(wx.StaticText(VectorRBDisplay,-1,
1862                    'Orientation reference atoms A-B-C: %s, %s, %s'%(atNames[rbRef[0]], \
1863                     atNames[rbRef[1]],atNames[rbRef[2]])),0)
1864            else:
1865                refAtmSizer.Add(wx.StaticText(VectorRBDisplay,-1,
1866                    'Orientation reference atoms A-B-C: '),0,wx.ALIGN_CENTER_VERTICAL)
1867                for i in range(3):
1868                    choices = [atNames[j] for j in refChoice[rbId][i]]
1869                    refSel = wx.ComboBox(VectorRBDisplay,-1,value='',
1870                        choices=choices,style=wx.CB_READONLY|wx.CB_DROPDOWN)
1871                    refSel.SetValue(atNames[rbRef[i]])
1872                    refSel.Bind(wx.EVT_COMBOBOX, OnRefSel)
1873                    Indx[refSel.GetId()] = i
1874                    refAtmSizer.Add(refSel,0,wx.ALIGN_CENTER_VERTICAL)
1875            return refAtmSizer
1876                       
1877        def rbVectMag(rbId,imag,rbData):
1878           
1879            def OnRBVectorMag(event):
1880                event.Skip()
1881                Obj = event.GetEventObject()
1882                rbId,imag = Indx[Obj.GetId()]
1883                try:
1884                    val = float(Obj.GetValue())
1885                    if val <= 0.:
1886                        raise ValueError
1887                    rbData['VectMag'][imag] = val
1888                except ValueError:
1889                    pass
1890                Obj.SetValue('%8.4f'%(val))
1891                wx.CallAfter(UpdateVectorRB,VectorRB.GetScrollPos(wx.VERTICAL))
1892                G2plt.PlotRigidBody(G2frame,'Vector',AtInfo,data['Vector'][rbId],plotDefaults)
1893               
1894            def OnRBVectorRef(event):
1895                Obj = event.GetEventObject()
1896                rbId,imag = Indx[Obj.GetId()]
1897                rbData['VectRef'][imag] = Obj.GetValue()
1898                       
1899            magSizer = wx.BoxSizer(wx.HORIZONTAL)
1900            magSizer.Add(wx.StaticText(VectorRBDisplay,-1,'Translation magnitude: '),
1901                0,wx.ALIGN_CENTER_VERTICAL)
1902            magValue = wx.TextCtrl(VectorRBDisplay,-1,'%8.4f'%(rbData['VectMag'][imag]))
1903            Indx[magValue.GetId()] = [rbId,imag]
1904            magValue.Bind(wx.EVT_TEXT_ENTER,OnRBVectorMag)
1905            magValue.Bind(wx.EVT_KILL_FOCUS,OnRBVectorMag)
1906            magSizer.Add(magValue,0,wx.ALIGN_CENTER_VERTICAL)
1907            magSizer.Add((5,0),)
1908            magref = wx.CheckBox(VectorRBDisplay,-1,label=' Refine?') 
1909            magref.SetValue(rbData['VectRef'][imag])
1910            magref.Bind(wx.EVT_CHECKBOX,OnRBVectorRef)
1911            Indx[magref.GetId()] = [rbId,imag]
1912            magSizer.Add(magref,0,wx.ALIGN_CENTER_VERTICAL)
1913            return magSizer
1914           
1915        def rbVectors(rbId,imag,mag,XYZ,rbData):
1916
1917            def TypeSelect(event):
1918                AtInfo = data['Vector']['AtInfo']
1919                r,c = event.GetRow(),event.GetCol()
1920                if vecGrid.GetColLabelValue(c) == 'Type':
1921                    PE = G2elemGUI.PickElement(G2frame,oneOnly=True)
1922                    if PE.ShowModal() == wx.ID_OK:
1923                        if PE.Elem != 'None':
1924                            El = PE.Elem.strip().lower().capitalize()
1925                            if El not in AtInfo:
1926                                Info = G2elem.GetAtomInfo(El)
1927                                AtInfo[El] = [Info['Drad'],Info['Color']]
1928                            rbData['rbTypes'][r] = El
1929                            vecGrid.SetCellValue(r,c,El)
1930                    PE.Destroy()
1931                wx.CallAfter(UpdateVectorRB,VectorRB.GetScrollPos(wx.VERTICAL))
1932
1933            def ChangeCell(event):
1934                r,c =  event.GetRow(),event.GetCol()
1935                if r >= 0 and (0 <= c < 3):
1936                    try:
1937                        val = float(vecGrid.GetCellValue(r,c))
1938                        rbData['rbVect'][imag][r][c] = val
1939                    except ValueError:
1940                        pass
1941                G2plt.PlotRigidBody(G2frame,'Vector',AtInfo,data['Vector'][rbId],plotDefaults)
1942                wx.CallAfter(UpdateVectorRB,VectorRB.GetScrollPos(wx.VERTICAL))
1943
1944            vecSizer = wx.BoxSizer()
1945            Types = 3*[wg.GRID_VALUE_FLOAT+':10,5',]+[wg.GRID_VALUE_STRING,]+3*[wg.GRID_VALUE_FLOAT+':10,5',]
1946            colLabels = ['Vector x','Vector y','Vector z','Type','Cart x','Cart y','Cart z']
1947            table = []
1948            rowLabels = []
1949            for ivec,xyz in enumerate(rbData['rbVect'][imag]):
1950                table.append(list(xyz)+[rbData['rbTypes'][ivec],]+list(XYZ[ivec]))
1951                rowLabels.append(str(ivec))
1952            vecTable = G2G.Table(table,rowLabels=rowLabels,colLabels=colLabels,types=Types)
1953            vecGrid = G2G.GSGrid(VectorRBDisplay)
1954            vecGrid.SetTable(vecTable, True)
1955            if 'phoenix' in wx.version():
1956                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGED, ChangeCell)
1957            else:
1958                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGE, ChangeCell)
1959            if not imag:
1960                vecGrid.Bind(wg.EVT_GRID_CELL_LEFT_DCLICK, TypeSelect)
1961            attr = wx.grid.GridCellAttr()
1962            attr.SetEditor(G2G.GridFractionEditor(vecGrid))
1963            for c in range(3):
1964                vecGrid.SetColAttr(c, attr)
1965            for row in range(vecTable.GetNumberRows()):
1966                if imag:
1967                    vecGrid.SetCellStyle(row,3,VERY_LIGHT_GREY,True)                   
1968                for col in [4,5,6]:
1969                    vecGrid.SetCellStyle(row,col,VERY_LIGHT_GREY,True)
1970#            vecGrid.SetScrollRate(0,0)
1971            vecGrid.AutoSizeColumns(False)
1972            vecSizer.Add(vecGrid)
1973            return vecSizer
1974       
1975        def FillRefChoice(rbId,rbData):
1976            choiceIds = [i for i in range(len(rbData['rbTypes']))]
1977           
1978            rbRef = rbData.get('rbRef',[-1,-1,-1,False])
1979            for i in range(3):
1980                choiceIds.remove(rbRef[i])
1981            refChoice[rbId] = [choiceIds[:],choiceIds[:],choiceIds[:]]
1982            for i in range(3):
1983                refChoice[rbId][i].append(rbRef[i])
1984                refChoice[rbId][i].sort()     
1985           
1986        if VectorRB.GetSizer(): VectorRB.GetSizer().Clear(True)
1987        VectorRBSizer = wx.BoxSizer(wx.VERTICAL)
1988        for rbId in data['RBIds']['Vector']:
1989            if rbId != 'AtInfo':
1990                rbData = data['Vector'][rbId]
1991                FillRefChoice(rbId,rbData)
1992                VectorRBSizer.Add(rbNameSizer(rbId,rbData),0)
1993                VectorRBSizer.Add(rbRefAtmSizer(rbId,rbData),0)
1994                XYZ = np.array([[0.,0.,0.] for Ty in rbData['rbTypes']])
1995                for imag,mag in enumerate(rbData['VectMag']):
1996                    XYZ += mag*rbData['rbVect'][imag]
1997                    VectorRBSizer.Add(rbVectMag(rbId,imag,rbData),0)
1998                    VectorRBSizer.Add(rbVectors(rbId,imag,mag,XYZ,rbData),0)
1999                VectorRBSizer.Add((5,5),0)
2000                data['Vector'][rbId]['rbXYZ'] = XYZ       
2001        VectorRBSizer.Layout()   
2002        VectorRBDisplay.SetSizer(VectorRBSizer,True)
2003        Size = VectorRBSizer.GetMinSize()
2004        Size[0] += 40
2005        Size[1] = max(Size[1],450) + 20
2006        VectorRBDisplay.SetSize(Size)
2007        VectorRB.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
2008        VectorRB.Scroll(0,Scroll)
2009       
2010    def UpdateResidueRB(rbId=0):
2011        '''Draw the contents of the Residue Rigid Body tab for Rigid Bodies tree entry
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            # start of rbNameSizer
2052            nameSizer = wx.BoxSizer(wx.HORIZONTAL)
2053            nameSizer.Add(wx.StaticText(ResidueRBDisplay,-1,'Residue name: '),
2054                0,wx.ALIGN_CENTER_VERTICAL)
2055            RBname = wx.TextCtrl(ResidueRBDisplay,-1,rbData['RBname'])
2056            Indx[RBname.GetId()] = rbId
2057            RBname.Bind(wx.EVT_TEXT_ENTER,OnRBName)
2058            RBname.Bind(wx.EVT_KILL_FOCUS,OnRBName)
2059            nameSizer.Add(RBname,0,wx.ALIGN_CENTER_VERTICAL)
2060            nameSizer.Add((5,0),)
2061            plotRB = wx.CheckBox(ResidueRBDisplay,-1,'Plot?')
2062            Indx[plotRB.GetId()] = rbId
2063            plotRB.Bind(wx.EVT_CHECKBOX,OnPlotRB)
2064            nameSizer.Add(plotRB,0,wx.ALIGN_CENTER_VERTICAL)
2065            nameSizer.Add((5,0),)
2066            if not rbData['useCount']:
2067                delRB = wx.CheckBox(ResidueRBDisplay,-1,'Delete?')
2068                Indx[delRB.GetId()] = rbId
2069                delRB.Bind(wx.EVT_CHECKBOX,OnDelRB)
2070                nameSizer.Add(delRB,0,wx.ALIGN_CENTER_VERTICAL)
2071                if 'H'  in rbData['rbTypes']:
2072                    stripH = wx.CheckBox(ResidueRBDisplay,-1,'Strip H-atoms?')
2073                    Indx[stripH.GetId()] = rbId
2074                    stripH.Bind(wx.EVT_CHECKBOX,OnStripH)
2075                    nameSizer.Add(stripH,0,wx.ALIGN_CENTER_VERTICAL)
2076            return nameSizer
2077           
2078        def rbResidues(rbId,rbData):
2079           
2080            def TypeSelect(event):
2081                AtInfo = data['Residue']['AtInfo']
2082                r,c = event.GetRow(),event.GetCol()
2083                if vecGrid.GetColLabelValue(c) == 'Type':
2084                    PE = G2elemGUI.PickElement(G2frame,oneOnly=True)
2085                    if PE.ShowModal() == wx.ID_OK:
2086                        if PE.Elem != 'None':
2087                            El = PE.Elem.strip().lower().capitalize()
2088                            if El not in AtInfo:
2089                                Info = G2elem.GetAtomInfo(El)
2090                                AtInfo[El] = [Info['Drad']['Color']]
2091                            rbData['rbTypes'][r] = El
2092                            vecGrid.SetCellValue(r,c,El)
2093                    PE.Destroy()
2094
2095            def ChangeCell(event):
2096                r,c =  event.GetRow(),event.GetCol()
2097                if r >= 0 and (0 <= c < 3):
2098                    try:
2099                        val = float(vecGrid.GetCellValue(r,c))
2100                        rbData['rbXYZ'][r][c] = val
2101                    except ValueError:
2102                        pass
2103                       
2104            def RowSelect(event):
2105                r,c =  event.GetRow(),event.GetCol()
2106                if c < 0:                   #only row clicks
2107                    for vecgrid in resList:
2108                        vecgrid.ClearSelection()
2109                    vecGrid.SelectRow(r,True)
2110
2111            def OnRefSel(event):
2112               
2113                Obj = event.GetEventObject()
2114                iref,res,jref = Indx[Obj.GetId()]
2115                sel = Obj.GetValue()
2116                ind = atNames.index(sel)
2117                if rbData['rbTypes'][ind] == 'H':
2118                    G2G.G2MessageBox(G2frame,'You should not select an H-atom for rigid body orientation')
2119                rbData['rbRef'][iref] = ind
2120                FillRefChoice(rbId,rbData)
2121                for i,ref in enumerate(RefObjs[jref]):
2122                    ref.SetItems([atNames[j] for j in refChoice[rbId][i]])
2123                    ref.SetValue(atNames[rbData['rbRef'][i]])                   
2124                rbXYZ = rbData['rbXYZ']
2125                if not iref:     #origin change
2126                    rbXYZ -= rbXYZ[ind]
2127                #TODO - transform all atom XYZ by axis choices
2128                Xxyz = rbXYZ[rbData['rbRef'][1]]
2129                X = Xxyz/np.sqrt(np.sum(Xxyz**2))
2130                Yxyz = rbXYZ[rbData['rbRef'][2]]
2131                Y = Yxyz/np.sqrt(np.sum(Yxyz**2))
2132                Mat = G2mth.getRBTransMat(X,Y)
2133                rbXYZ = np.inner(Mat,rbXYZ).T
2134                rbData['rbXYZ'] = rbXYZ
2135                res.ClearSelection()
2136                resTable = res.GetTable()
2137                for r in range(res.GetNumberRows()):
2138                    row = resTable.GetRowValues(r)
2139                    row[2:4] = rbXYZ[r]
2140                    resTable.SetRowValues(r,row)
2141                res.ForceRefresh()
2142                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
2143               
2144            Types = 2*[wg.GRID_VALUE_STRING,]+3*[wg.GRID_VALUE_FLOAT+':10,5',]
2145            colLabels = ['Name','Type','Cart x','Cart y','Cart z']
2146            table = []
2147            rowLabels = []
2148            for ivec,xyz in enumerate(rbData['rbXYZ']):
2149                table.append([rbData['atNames'][ivec],]+[rbData['rbTypes'][ivec],]+list(xyz))
2150                rowLabels.append(str(ivec))
2151            vecTable = G2G.Table(table,rowLabels=rowLabels,colLabels=colLabels,types=Types)
2152            vecGrid = G2G.GSGrid(ResidueRBDisplay)
2153            Indx[vecGrid.GetId()] = rbId
2154            resList.append(vecGrid)
2155            vecGrid.SetTable(vecTable, True)
2156            if 'phoenix' in wx.version():
2157                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGED, ChangeCell)
2158            else:
2159                vecGrid.Bind(wg.EVT_GRID_CELL_CHANGE, ChangeCell)
2160            vecGrid.Bind(wg.EVT_GRID_CELL_LEFT_DCLICK, TypeSelect)
2161            vecGrid.Bind(wg.EVT_GRID_LABEL_LEFT_CLICK, RowSelect)
2162            attr = wx.grid.GridCellAttr()
2163            attr.SetEditor(G2G.GridFractionEditor(vecGrid))
2164            for c in range(3):
2165                vecGrid.SetColAttr(c, attr)
2166            for row in range(vecTable.GetNumberRows()):
2167                for col in range(5):
2168                    vecGrid.SetCellStyle(row,col,VERY_LIGHT_GREY,True)
2169            vecGrid.AutoSizeColumns(False)
2170            vecSizer = wx.BoxSizer()
2171            vecSizer.Add(vecGrid)
2172           
2173            refAtmSizer = wx.BoxSizer(wx.HORIZONTAL)
2174            atNames = rbData['atNames']
2175            rbRef = rbData['rbRef']
2176            if rbData['rbRef'][3] or rbData['useCount']:
2177                refAtmSizer.Add(wx.StaticText(ResidueRBDisplay,-1,
2178                    'Orientation reference non-H atoms A-B-C: %s, %s, %s'%(atNames[rbRef[0]], \
2179                     atNames[rbRef[1]],atNames[rbRef[2]])),0)
2180            else:
2181                refAtmSizer.Add(wx.StaticText(ResidueRBDisplay,-1,
2182                    'Orientation reference non-H atoms A-B-C: '),0,wx.ALIGN_CENTER_VERTICAL)
2183                refObj = [0,0,0]
2184                for i in range(3):
2185                    choices = [atNames[j] for j in refChoice[rbId][i]]
2186                    refSel = wx.ComboBox(ResidueRBDisplay,-1,value='',
2187                        choices=choices,style=wx.CB_READONLY|wx.CB_DROPDOWN)
2188                    refSel.SetValue(atNames[rbRef[i]])
2189                    refSel.Bind(wx.EVT_COMBOBOX, OnRefSel)
2190                    Indx[refSel.GetId()] = [i,vecGrid,len(RefObjs)]
2191                    refObj[i] = refSel
2192                    refAtmSizer.Add(refSel,0,wx.ALIGN_CENTER_VERTICAL)
2193                RefObjs.append(refObj)
2194           
2195            mainSizer = wx.BoxSizer(wx.VERTICAL)
2196            mainSizer.Add(refAtmSizer)
2197            mainSizer.Add(vecSizer)
2198            return mainSizer
2199           
2200        def SeqSizer(angSlide,rbId,iSeq,Seq,atNames):
2201           
2202            def ChangeAngle(event):
2203                event.Skip()
2204                Obj = event.GetEventObject()
2205                rbId,Seq = Indx[Obj.GetId()][:2]
2206                val = Seq[2]
2207                try:
2208                    val = float(Obj.GetValue())
2209                    Seq[2] = val
2210                except ValueError:
2211                    pass
2212                Obj.SetValue('%8.2f'%(val))
2213                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,data['Residue'][rbId],plotDefaults)
2214               
2215            def OnRadBtn(event):
2216                Obj = event.GetEventObject()
2217                Seq,iSeq,angId = Indx[Obj.GetId()]
2218                data['Residue'][rbId]['SelSeq'] = [iSeq,angId]
2219                angSlide.SetValue(int(100*Seq[2]))
2220               
2221            def OnDelBtn(event):
2222                Obj = event.GetEventObject()
2223                rbId,Seq = Indx[Obj.GetId()]
2224                data['Residue'][rbId]['rbSeq'].remove(Seq)       
2225                wx.CallAfter(UpdateResidueRB,rbId)
2226           
2227            seqSizer = wx.FlexGridSizer(0,5,2,2)
2228            seqSizer.AddGrowableCol(3,0)
2229            iBeg,iFin,angle,iMove = Seq
2230            ang = wx.TextCtrl(ResidueRBDisplay,-1,'%8.2f'%(angle),size=(50,20))
2231            if not iSeq:
2232                radBt = wx.RadioButton(ResidueRBDisplay,-1,'',style=wx.RB_GROUP)
2233                data['Residue'][rbId]['SelSeq'] = [iSeq,ang.GetId()]
2234            else:
2235                radBt = wx.RadioButton(ResidueRBDisplay,-1,'')
2236            radBt.Bind(wx.EVT_RADIOBUTTON,OnRadBtn)                   
2237            seqSizer.Add(radBt)
2238            delBt = wx.RadioButton(ResidueRBDisplay,-1,'')
2239            delBt.Bind(wx.EVT_RADIOBUTTON,OnDelBtn)
2240            seqSizer.Add(delBt)
2241            bond = wx.TextCtrl(ResidueRBDisplay,-1,'%s %s'%(atNames[iBeg],atNames[iFin]),size=(50,20))
2242            seqSizer.Add(bond,0,wx.ALIGN_CENTER_VERTICAL)
2243            Indx[radBt.GetId()] = [Seq,iSeq,ang.GetId()]
2244            Indx[delBt.GetId()] = [rbId,Seq]
2245            Indx[ang.GetId()] = [rbId,Seq,ang]
2246            ang.Bind(wx.EVT_TEXT_ENTER,ChangeAngle)
2247            ang.Bind(wx.EVT_KILL_FOCUS,ChangeAngle)
2248            seqSizer.Add(ang,0,wx.ALIGN_CENTER_VERTICAL)
2249            atms = ''
2250            for i in iMove:   
2251                atms += ' %s,'%(atNames[i])
2252            moves = wx.TextCtrl(ResidueRBDisplay,-1,atms[:-1],size=(200,20))
2253            seqSizer.Add(moves,1,wx.ALIGN_CENTER_VERTICAL|wx.EXPAND|wx.RIGHT)
2254            return seqSizer
2255           
2256        def SlideSizer():
2257           
2258            def OnSlider(event):
2259                Obj = event.GetEventObject()
2260                rbData = Indx[Obj.GetId()]
2261                iSeq,angId = rbData['SelSeq']
2262                val = float(Obj.GetValue())/100.
2263                rbData['rbSeq'][iSeq][2] = val
2264                Indx[angId][2].SetValue('%8.2f'%(val))
2265                G2plt.PlotRigidBody(G2frame,'Residue',AtInfo,rbData,plotDefaults)
2266           
2267            slideSizer = wx.BoxSizer(wx.HORIZONTAL)
2268            slideSizer.Add(wx.StaticText(ResidueRBDisplay,-1,'Selected torsion angle:'),0)
2269            iSeq,angId = rbData['SelSeq']
2270            angSlide = wx.Slider(ResidueRBDisplay,-1,
2271                int(100*rbData['rbSeq'][iSeq][2]),0,36000,size=(200,20),
2272                style=wx.SL_HORIZONTAL)
2273            angSlide.Bind(wx.EVT_SLIDER, OnSlider)
2274            Indx[angSlide.GetId()] = rbData
2275            slideSizer.Add(angSlide,0)           
2276            return slideSizer,angSlide
2277           
2278        def FillRefChoice(rbId,rbData):
2279            choiceIds = [i for i in range(len(rbData['atNames']))]
2280            for seq in rbData['rbSeq']:
2281                for i in seq[3]:
2282                    try:
2283                        choiceIds.remove(i)
2284                    except ValueError:
2285                        pass
2286            rbRef = rbData['rbRef']
2287            for i in range(3):
2288                try:
2289                    choiceIds.remove(rbRef[i])
2290                except ValueError:
2291                    pass
2292            refChoice[rbId] = [choiceIds[:],choiceIds[:],choiceIds[:]]
2293            for i in range(3):
2294                refChoice[rbId][i].append(rbRef[i])
2295                refChoice[rbId][i].sort()
2296               
2297        def OnSelect(event):
2298            rbname = rbchoice[select.GetSelection()]
2299            rbId = RBnames[rbname]
2300            wx.CallLater(100,UpdateResidueRB,rbId)
2301           
2302        #----- beginning of UpdateResidueRB -----------------------------------------------
2303        AtInfo = data['Residue']['AtInfo']
2304        refChoice = {}
2305        RefObjs = []
2306
2307        GS = ResidueRBDisplay.GetSizer()
2308        if GS: 
2309            try:        #get around a c++ error in wx 4.0; doing is again seems to be OK
2310                GS.Clear(True)
2311            except:
2312                GS.Clear(True)
2313       
2314        RBnames = {}
2315        for rbid in data['RBIds']['Residue']:
2316            RBnames.update({data['Residue'][rbid]['RBname']:rbid,})
2317        if not RBnames:
2318            return
2319        rbchoice = list(RBnames.keys())
2320        if GS: 
2321            ResidueRBSizer = GS
2322        else:
2323            ResidueRBSizer = wx.BoxSizer(wx.VERTICAL)
2324        if len(RBnames) > 1:
2325            selSizer = wx.BoxSizer(wx.HORIZONTAL)
2326            selSizer.Add(wx.StaticText(ResidueRBDisplay,label=' Select residue to view:'),0)
2327            rbchoice.sort()
2328            select = wx.ComboBox(ResidueRBDisplay,choices=rbchoice)
2329            select.Bind(wx.EVT_COMBOBOX,OnSelect)
2330            selSizer.Add(select,0)
2331            ResidueRBSizer.Add(selSizer,0)
2332        if not rbId:
2333            rbId = RBnames[rbchoice[0]]
2334        rbData = data['Residue'][rbId]
2335        FillRefChoice(rbId,rbData)
2336        ResidueRBSizer.Add(rbNameSizer(rbId,rbData),0)
2337        ResidueRBSizer.Add(rbResidues(rbId,rbData),0)
2338        ResidueRBSizer.Add((5,5),0)
2339        if rbData['rbSeq']:
2340            slideSizer,angSlide = SlideSizer()
2341        if len(rbData['rbSeq']):
2342            ResidueRBSizer.Add(wx.StaticText(ResidueRBDisplay,-1,
2343                'Sel  Del  Bond             Angle      Riding atoms'),
2344                0,wx.ALIGN_CENTER_VERTICAL)                       
2345        for iSeq,Seq in enumerate(rbData['rbSeq']):
2346            ResidueRBSizer.Add(SeqSizer(angSlide,rbId,iSeq,Seq,rbData['atNames']))
2347        if rbData['rbSeq']:
2348            ResidueRBSizer.Add(slideSizer,)
2349
2350        ResidueRBSizer.Add((5,25),)
2351        ResidueRBSizer.Layout()   
2352        ResidueRBDisplay.SetSizer(ResidueRBSizer,True)
2353        ResidueRBDisplay.SetAutoLayout(True)
2354        Size = ResidueRBSizer.GetMinSize()
2355        ResidueRBDisplay.SetSize(Size)
2356       
2357        Size[0] += 40
2358        Size[1] = max(Size[1],450) + 20
2359        ResidueRB.SetSize(Size)
2360        ResidueRB.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)       
2361       
2362        ResidueRBDisplay.Show()
2363       
2364    def SetStatusLine(text):
2365        G2frame.GetStatusBar().SetStatusText(text,1)                                     
2366
2367    G2gd.SetDataMenuBar(G2frame,G2frame.dataWindow.RigidBodyMenu)
2368    SetStatusLine('')
2369    UpdateVectorRB()
2370    G2frame.rbBook.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, OnPageChanged)
2371    wx.CallAfter(OnPageChanged,None)
Note: See TracBrowser for help on using the repository browser.