Changeset 819


Ignore:
Timestamp:
Dec 21, 2012 4:50:13 PM (9 years ago)
Author:
vondreele
Message:

restraint fixes

Location:
trunk
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIrestrGUI.py

    r818 r819  
    3939    restrData = data[phaseName]
    4040    if 'Bond' not in restrData:
    41         restrData['Bond'] = {'wtFactor':1.0,'Range':0.9,'Bonds':[],'Use':True}
     41        restrData['Bond'] = {'wtFactor':1.0,'Range':1.1,'Bonds':[],'Use':True}
    4242    if 'Angle' not in restrData:
    4343        restrData['Angle'] = {'wtFactor':1.0,'Range':0.85,'Angles':[],'Use':True}
    4444    if 'Plane' not in restrData:
    45         restrData['Plane'] = {'wtFactor':1.0,'Range':0.9,'Planes':[],'Use':True}
     45        restrData['Plane'] = {'wtFactor':1.0,'Planes':[],'Use':True}
    4646    if 'Chiral' not in restrData:
    47         restrData['Chiral'] = {'wtFactor':1.0,'Range':0.9,'Volumes':[],'Use':True}
     47        restrData['Chiral'] = {'wtFactor':1.0,'Volumes':[],'Use':True}
    4848    if 'Torsion' not in restrData:
    49         restrData['Torsion'] = {'wtFactor':1.0,'Range':0.9,'Coeff':{},'Torsions':[],'Use':True}
     49        restrData['Torsion'] = {'wtFactor':1.0,'Coeff':{},'Torsions':[],'Use':True}
    5050    if 'Rama' not in restrData:
    51         restrData['Rama'] = {'wtFactor':1.0,'Range':0.9,'Coeff':{},'Ramas':[],'Use':True}
     51        restrData['Rama'] = {'wtFactor':1.0,'Coeff':{},'Ramas':[],'Use':True}
    5252    General = phasedata['General']
    5353    Cell = General['Cell'][1:7]          #skip flag & volume   
     
    168168        elif 'Plane' in G2frame.dataDisplay.GetPageText(page):
    169169            AddPlaneRestraint(restrData['Plane'])
    170         elif 'Chiral' in G2frame.dataDisplay.GetPageText(page):
    171             AddChiralRestraint(restrData['Chiral'])
    172 #        elif 'Torsion' in G2frame.dataDisplay.GetPageText(page):
    173 #            AddTorsionRestraint(restrData['Torsion'])
    174 #        elif 'Rama' in G2frame.dataDisplay.GetPageText(page):
    175 #            AddRamaRestraint(restrData['Rama'])
    176170           
    177171    def OnAddAARestraint(event):
     
    191185           
    192186    def AddBondRestraint(bondRestData):
    193         Radii = dict(zip(General['AtomTypes'],General['BondRadii']))
    194187        Lists = {'origin':[],'target':[]}
    195188        for listName in ['origin','target']:
     
    206199                    else:
    207200                        Lists[listName].append([Ids[x],Types[x],Coords[x],])
     201        bond = 1.54
     202        dlg = G2phG.SingleFloatDialog(G2frame,'Distance','Enter restraint distance for bond',bond,[1.,4.],'%.4f')
     203        if dlg.ShowModal() == wx.ID_OK:
     204            bond = dlg.GetValue()
     205        dlg.Destroy()
    208206        Factor = bondRestData['Range']
    209207        indices = (-1,0,1)
     
    223221                    else:
    224222                        result = G2spc.GenAtom(Tcoord,SGData,False,Move=False)
    225                     BsumR = (Radii[Otype]+Radii[Ttype])*Factor
    226223                    for Txyz,Top,Tunit in result:
    227224                        Dx = (Txyz-np.array(Ocoord))+Units
    228225                        dx = np.inner(Amat,Dx)
    229                         dist = ma.masked_less(np.sqrt(np.sum(dx**2,axis=0)),0.5)
    230                         IndB = ma.nonzero(ma.masked_greater(dist-BsumR,0.))
     226                        dist = ma.masked_less(np.sqrt(np.sum(dx**2,axis=0)),bond/Factor)
     227                        IndB = ma.nonzero(ma.masked_greater(dist,bond*Factor))
    231228                        if np.any(IndB):
    232229                            for indb in IndB:
     
    237234                                    else:
    238235                                        Topstr = str(Top)
    239                                     newBond = [[Oid,Tid],['1',Topstr],1.54,0.01]
     236                                    newBond = [[Oid,Tid],['1',Topstr],bond,0.01]
    240237                                    if newBond not in bondRestData['Bonds']:
    241238                                        bondRestData['Bonds'].append(newBond)
     
    278275                    for tId,tCoord in zip(tIds,tCoords)[i:]:
    279276                        obsd = np.sqrt(np.sum(np.inner(Amat,tCoord-oCoord)**2))
    280                         if dist*Factor < obsd < dist/Factor:
     277                        if dist/Factor < obsd < dist*Factor:
    281278                            newBond = [[oId,tId],['1','1'],dist,esd]
    282279                            if newBond not in bondRestData['Bonds']:
     
    302299                    origAtoms.append([Ids[x],Types[x],Coords[x]])
    303300        targAtoms = [[Ids[x+iBeg],Types[x+iBeg],Coords[x+iBeg]] for x in range(len(Names[iBeg:]))]
     301        value = 109.54
     302        dlg = G2phG.SingleFloatDialog(G2frame,'Angle','Enter restraint angle ',value,[30.,180.],'%.2f')
     303        if dlg.ShowModal() == wx.ID_OK:
     304            value = dlg.GetValue()
     305        dlg.Destroy()
    304306
    305307        Factor = angleRestData['Range']
     
    340342                        ops = [vecta[4],vecta[1],vectb[4]]
    341343                        XYZ = np.array([vecta[5],vecta[2],vectb[5]])
    342                         angle = G2mth.getRestAngle(XYZ,Amat)
     344                        angle = [ids,ops,value,1.0]
    343345                        if angle not in angleRestData['Angles']:
    344                             angleRestData['Angles'].append([ids,ops,109.5,1.0])
     346                            angleRestData['Angles'].append(angle)
    345347        UpdateAngleRestr(angleRestData)               
    346348
     
    411413                    restrData['Planes'].append(plane)
    412414            else:
    413                 print '**** ERROR - not enough atoms fro a plane restraint - try again ****'
     415                print '**** ERROR - not enough atoms for a plane restraint - try again ****'
    414416        UpdatePlaneRestr(restrData)               
    415417
     
    466468        UpdatePlaneRestr(planeRestData)               
    467469
    468     def AddChiralRestraint():
    469         print 'Chiral restraint'
    470        
    471470    def AddAAChiralRestraint(chiralRestData):
    472471        macro = getMacroFile('chiral')
     
    517516        return Chains
    518517       
    519 #       
    520518    def AddAATorsionRestraint(torsionRestData):
    521519        macro = getMacroFile('torsion')
     
    665663               
    666664    def WtBox(wind,restData):
    667         if 'Range' not in restData: restData['Range'] = 0.9     #patch
     665        if 'Range' not in restData: restData['Range'] = 1.1     #patch
    668666       
    669667        def OnWtFactor(event):
     
    696694        useData.SetValue(restData['Use'])       
    697695        wtBox.Add(useData,0,wx.ALIGN_CENTER_VERTICAL)
    698         if 'Bond' in restData or 'Angle' in restData:
     696        if 'Bonds' in restData or 'Angles' in restData:
    699697            wtBox.Add(wx.StaticText(wind,-1,'Search range:'),0,wx.ALIGN_CENTER_VERTICAL)
    700698            sRange = wx.TextCtrl(wind,-1,value='%.2f'%(restData['Range']),style=wx.TE_PROCESS_ENTER)
     
    738736                return
    739737            Bonds.ClearSelection()
    740             val = bondList[rows[0]][3]
     738            val = bondList[rows[0]][2]
    741739            dlg = G2phG.SingleFloatDialog(G2frame,'New value','Enter new value for bond',val,[0.,5.],'%.4f')
    742740            if dlg.ShowModal() == wx.ID_OK:
    743741                parm = dlg.GetValue()
    744742                for r in rows:
    745                     bondList[r][3] = parm
     743                    bondRestData['Bonds'][r][2] = parm
    746744            dlg.Destroy()
    747745            UpdateBondRestr(bondRestData)               
     
    752750                return
    753751            Bonds.ClearSelection()
    754             val = bondList[rows[0]][4]
     752            val = bondList[rows[0]][3]
    755753            dlg = G2phG.SingleFloatDialog(G2frame,'New value','Enter new esd for bond',val,[0.,1.],'%.4f')
    756754            if dlg.ShowModal() == wx.ID_OK:
    757755                parm = dlg.GetValue()
    758756                for r in rows:
    759                     bondList[r][4] = parm
     757                    bondRestData['Bonds'][r][3] = parm
    760758            dlg.Destroy()
    761759            UpdateBondRestr(bondRestData)               
     
    798796            else:
    799797                colLabels = ['A+SymOp - B+SymOp','calc','obs','esd']
    800                 for i,[indx,ops,obs,esd] in enumerate(bondList):
    801                     names = G2mth.GetAtomItemsById(Atoms,AtLookUp,indx,ct-1)
    802                     XYZ = np.array(G2mth.GetAtomItemsById(Atoms,AtLookUp,indx,cx,3))
    803                     XYZ = G2mth.getSyXYZ(XYZ,ops,SGData)
    804                     calc = G2mth.getRestDist(XYZ,Amat)
    805                     table.append([names[0]+'+('+ops[0]+') - '+names[1]+'+('+ops[1]+')',calc,obs,esd])
    806                     rowLabels.append(str(i))
     798                try:
     799                    for i,[indx,ops,obs,esd] in enumerate(bondList):
     800                        names = G2mth.GetAtomItemsById(Atoms,AtLookUp,indx,ct-1)
     801                        XYZ = np.array(G2mth.GetAtomItemsById(Atoms,AtLookUp,indx,cx,3))
     802                        XYZ = G2mth.getSyXYZ(XYZ,ops,SGData)
     803                        calc = G2mth.getRestDist(XYZ,Amat)
     804                        table.append([names[0]+'+('+ops[0]+') - '+names[1]+'+('+ops[1]+')',calc,obs,esd])
     805                        rowLabels.append(str(i))
     806                except ValueError:              #patch - old restraints!
     807                    bondList = []
    807808            bondTable = G2gd.Table(table,rowLabels=rowLabels,colLabels=colLabels,types=Types)
    808809            Bonds = G2gd.GSGrid(BondRestr)
     
    836837                return
    837838            Angles.ClearSelection()
    838             val = angleList[rows[0]][3]
     839            val = angleList[rows[0]][2]
    839840            dlg = G2phG.SingleFloatDialog(G2frame,'New value','Enter new value for angle',val,[0.,360.],'%.2f')
    840841            if dlg.ShowModal() == wx.ID_OK:
    841842                parm = dlg.GetValue()
    842843                for r in rows:
    843                     angleList[r][3] = parm
     844                    angleRestData['Angles'][r][2] = parm
    844845            dlg.Destroy()
    845846            UpdateAngleRestr(angleRestData)               
     
    850851                return
    851852            Angles.ClearSelection()
    852             val = angleList[rows[0]][4]
     853            val = angleList[rows[0]][3]
    853854            dlg = G2phG.SingleFloatDialog(G2frame,'New value','Enter new esd for angle',val,[0.,5.],'%.2f')
    854855            if dlg.ShowModal() == wx.ID_OK:
    855856                parm = dlg.GetValue()
    856857                for r in rows:
    857                     angleList[r][4] = parm
     858                    angleRestData['Angles'][r][3] = parm
    858859            dlg.Destroy()
    859860            UpdateAngleRestr(angleRestData)               
     
    925926        AngleRestr.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
    926927        G2frame.dataFrame.SetSize(Size)
    927 #        G2frame.dataFrame.setSizePosLeft(Size)
    928928   
    929929    def UpdatePlaneRestr(planeRestData):
     
    939939                return
    940940            Planes.ClearSelection()
    941             val = planeList[rows[0]][4]
     941            val = planeList[rows[0]][3]
    942942            dlg = G2phG.SingleFloatDialog(G2frame,'New value','Enter new esd for plane',val,[0.,5.],'%.2f')
    943943            if dlg.ShowModal() == wx.ID_OK:
    944944                parm = dlg.GetValue()
    945945                for r in rows:
    946                     planeList[r][4] = parm
     946                    planeRestData['Planes'][r][3] = parm
    947947            dlg.Destroy()
    948948            UpdatePlaneRestr(planeRestData)               
     
    10381038                return
    10391039            Volumes.ClearSelection()
    1040             val = volumeList[rows[0]][3]
     1040            val = volumeList[rows[0]][2]
    10411041            dlg = G2phG.SingleFloatDialog(G2frame,'New value','Enter new value for chiral volume',val,[0.,360.],'%.2f')
    10421042            if dlg.ShowModal() == wx.ID_OK:
    10431043                parm = dlg.GetValue()
    10441044                for r in rows:
    1045                     volumeList[r][3] = parm
     1045                    chiralRestData['Volumes'][r][2] = parm
    10461046            dlg.Destroy()
    10471047            UpdateChiralRestr(chiralRestData)               
     
    10521052                return
    10531053            Volumes.ClearSelection()
    1054             val = volumeList[rows[0]][4]
     1054            val = volumeList[rows[0]][3]
    10551055            dlg = G2phG.SingleFloatDialog(G2frame,'New value','Enter new esd for chiral volume',val,[0.,5.],'%.2f')
    10561056            if dlg.ShowModal() == wx.ID_OK:
    10571057                parm = dlg.GetValue()
    10581058                for r in rows:
    1059                     volumeList[r][4] = parm
     1059                    chiralRestData['Volumes'][r][3] = parm
    10601060            dlg.Destroy()
    10611061            UpdateChiralRestr(chiralRestData)               
     
    11171117        ChiralRestr.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
    11181118        G2frame.dataFrame.SetSize(Size)
    1119 #        G2frame.dataFrame.setSizePosLeft(Size)
    11201119   
    11211120    def UpdateTorsionRestr(torsionRestData):
     
    11411140                parm = dlg.GetValue()
    11421141                for r in rows:
    1143                     volumeList[r][4] = parm
     1142                    torsionRestData['Torsions'][r][4] = parm
    11441143            dlg.Destroy()
    11451144            UpdateTorsionRestr(torsionRestData)               
     
    12111210        TorsionRestr.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
    12121211        G2frame.dataFrame.SetSize(Size)
    1213 #        G2frame.dataFrame.setSizePosLeft(Size)
    12141212
    12151213    def UpdateRamaRestr(ramaRestData):
     
    12301228                return
    12311229            Ramas.ClearSelection()
    1232             val = ramaList[rows[0]][5]
     1230            val = ramaList[rows[0]][4]
    12331231            dlg = G2phG.SingleFloatDialog(G2frame,'New value','Enter new esd for energy',val,[0.,5.],'%.2f')
    12341232            if dlg.ShowModal() == wx.ID_OK:
    12351233                parm = dlg.GetValue()
    12361234                for r in rows:
    1237                     ramaList[r][4] = parm
     1235                    ramaRestData['Ramas'][r][4] = parm
    12381236            dlg.Destroy()
    12391237            UpdateRamaRestr(ramaRestData)               
     
    13311329        elif text == 'Chiral restraints':
    13321330            G2gd.SetDataMenuBar(G2frame,G2frame.dataFrame.RestraintMenu)
    1333             G2frame.dataFrame.RestraintEdit.Enable(G2gd.wxID_RESTRAINTADD,True)
     1331            G2frame.dataFrame.RestraintEdit.Enable(G2gd.wxID_RESTRAINTADD,False)
    13341332            G2frame.dataFrame.RestraintEdit.Enable(G2gd.wxID_RESRCHANGEVAL,True)
    13351333            chiralRestData = restrData['Chiral']
  • trunk/GSASIIstruct.py

    r811 r819  
    603603        if phaseRest:
    604604            print >>pFile,'\n Restraints:'
    605             names = ['Bonds','Angles','Planes','Volumes']
    606             for i,name in enumerate(['Bond','Angle','Plane','Chiral']):
     605            names = ['Bonds','Angles','Planes','Volumes','Torsions','Ramas']
     606            for i,name in enumerate(['Bond','Angle','Plane','Chiral','Torsion','Rama']):
    607607                itemRest = phaseRest[name]
    608608                if itemRest[names[i]]:
    609                     print >>pFile,'\n  %30s %10.3f'%(name+' restraint weight factor',itemRest['wtFactor'])
     609                    print >>pFile,'\n  %30s %10.3f Use: %s'%(name+' restraint weight factor',itemRest['wtFactor'],str(itemRest['Use']))
    610610                    print >>pFile,'  atoms(symOp), calc, obs, sig: '
    611611                    for item in phaseRest[name][names[i]]:
    612                         text = '   '
    613                         for a,at in enumerate(item[0]):
    614                             text += at+'+('+item[1][a]+') '
    615                             if (a+1)%5 == 0:
    616                                 text += '\n   '
    617                         print >>pFile,text,' %.3f %.3f %.3f'%(item[3],item[4],item[5])
     612                        print item
     613#                        text = '   '
     614#                        for a,at in enumerate(item[0]):
     615#                            text += at+'+('+item[1][a]+') '
     616#                            if (a+1)%5 == 0:
     617#                                text += '\n   '
     618#                        print >>pFile,text,' %.3f %.3f %.3f'%(item[3],item[4],item[5])
    618619       
    619620    if Print:print  >>pFile,' Phases:'
Note: See TracChangeset for help on using the changeset viewer.