Changeset 1951 for trunk/GSASIIphsGUI.py


Ignore:
Timestamp:
Aug 7, 2015 11:54:12 AM (8 years ago)
Author:
vondreele
Message:

implement drawing of incommensurate structures - live variation wrt tau
required careful attention to siteSym text & some reorganization of modules
more work on special pos code

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIphsGUI.py

    r1947 r1951  
    7878    mainFrame.SendSizeEvent()
    7979   
     80#def FindBondsToo():                         #works but slow for large structures - keep as reference
     81#    cx,ct,cs,ci = data['Drawing']['atomPtrs']
     82#    atomData = data['Drawing']['Atoms']
     83#    generalData = data['General']
     84#    Amat,Bmat = G2lat.cell2AB(generalData['Cell'][1:7])
     85#    radii = generalData['BondRadii']
     86#    atomTypes = generalData['AtomTypes']
     87#    try:
     88#        indH = atomTypes.index('H')
     89#        radii[indH] = 0.5
     90#    except:
     91#        pass           
     92#    for atom in atomData:
     93#        atom[-1] = []
     94#    Atoms = []
     95#    for i,atom in enumerate(atomData):
     96#        Atoms.append([i,np.array(atom[cx:cx+3]),atom[cs],radii[atomTypes.index(atom[ct])]])
     97#    for atomA in Atoms:
     98#        if atomA[2] in ['lines','sticks','ellipsoids','balls & sticks','polyhedra']:
     99#            for atomB in Atoms:                   
     100#                Dx = atomB[1]-atomA[1]
     101#                DX = np.inner(Amat,Dx)
     102#                dist = np.sqrt(np.sum(DX**2))
     103#                sumR = atomA[3]+atomB[3]
     104#                if 0.5 < dist <= 0.85*sumR:
     105#                    i = atomA[0]
     106#                    if atomA[2] == 'polyhedra':
     107#                        atomData[i][-1].append(DX)
     108#                    elif atomB[1] != 'polyhedra':
     109#                        j = atomB[0]
     110#                        atomData[i][-1].append(Dx*atomA[3]/sumR)
     111#                        atomData[j][-1].append(-Dx*atomB[3]/sumR)
     112               
     113def FindBondsDraw(data):   
     114    '''uses numpy & masks - very fast even for proteins!
     115    '''
     116    import numpy.ma as ma
     117    cx,ct,cs,ci = data['Drawing']['atomPtrs']
     118    hydro = data['Drawing']['showHydrogen']
     119    atomData = data['Drawing']['Atoms']
     120    generalData = data['General']
     121    Amat,Bmat = G2lat.cell2AB(generalData['Cell'][1:7])
     122    radii = generalData['BondRadii']
     123    atomTypes = generalData['AtomTypes']
     124    try:
     125        indH = atomTypes.index('H')
     126        radii[indH] = 0.5
     127    except:
     128        pass           
     129    for atom in atomData:
     130        atom[-2] = []               #clear out old bonds/polyhedra
     131        atom[-1] = []
     132    Indx = range(len(atomData))
     133    Atoms = []
     134    Styles = []
     135    Radii = []
     136    for atom in atomData:
     137        Atoms.append(np.array(atom[cx:cx+3]))
     138        Styles.append(atom[cs])
     139        try:
     140            if not hydro and atom[ct] == 'H':
     141                Radii.append(0.0)
     142            else:
     143                Radii.append(radii[atomTypes.index(atom[ct])])
     144        except ValueError:          #changed atom type!
     145            Radii.append(0.20)
     146    Atoms = np.array(Atoms)
     147    Radii = np.array(Radii)
     148    IASR = zip(Indx,Atoms,Styles,Radii)
     149    for atomA in IASR:
     150        if atomA[2] in ['lines','sticks','ellipsoids','balls & sticks','polyhedra']:
     151            Dx = Atoms-atomA[1]
     152            dist = ma.masked_less(np.sqrt(np.sum(np.inner(Amat,Dx)**2,axis=0)),0.5) #gets rid of G2frame & disorder "bonds" < 0.5A
     153            sumR = atomA[3]+Radii
     154            IndB = ma.nonzero(ma.masked_greater(dist-data['Drawing']['radiusFactor']*sumR,0.))                 #get indices of bonded atoms
     155            i = atomA[0]
     156            for j in IndB[0]:
     157                if Styles[i] == 'polyhedra':
     158                    atomData[i][-2].append(np.inner(Amat,Dx[j]))
     159                elif Styles[j] != 'polyhedra' and j > i:
     160                    atomData[i][-2].append(Dx[j]*Radii[i]/sumR[j])
     161                    atomData[j][-2].append(-Dx[j]*Radii[j]/sumR[j])
     162            if Styles[i] == 'polyhedra':
     163                Bonds = atomData[i][-2]
     164                Faces = []
     165                if len(Bonds) > 2:
     166                    FaceGen = G2lat.uniqueCombinations(Bonds,3)     #N.B. this is a generator
     167                    for face in FaceGen:
     168                        vol = nl.det(face)
     169                        if abs(vol) > 1. or len(Bonds) == 3:
     170                            if vol < 0.:
     171                                face = [face[0],face[2],face[1]]
     172                            face = np.array(face)
     173                            if not np.array([np.array(nl.det(face-bond))+0.0001 < 0 for bond in Bonds]).any():
     174                                norm = np.cross(face[1]-face[0],face[2]-face[0])
     175                                norm /= np.sqrt(np.sum(norm**2))
     176                                Faces.append([face,norm])
     177                    atomData[i][-1] = Faces
     178                       
    80179def UpdatePhaseData(G2frame,Item,data,oldPage):
    81180    '''Create the data display window contents when a phase is clicked on
     
    98197    def GetReflData(G2frame,phaseName,reflNames):
    99198        ReflData = {'RefList':[],'Type':''}
     199        if '' in reflNames:
     200            return None
    100201        for reflName in reflNames:
    101202            if 'PWDR' in reflName:
     
    12091310                                atomData[r][c] = parms
    12101311                        if 'Atoms' in data['Drawing']:
    1211                             DrawAtomsReplaceByID(data['Drawing'],atomData[r],ID)
     1312                            DrawAtomsReplaceByID(data['Drawing'],ui+6,atomData[r],ID)
    12121313                    wx.CallAfter(Paint)
    12131314                   
     
    12771378                    atomData[r][c] = atomData[r][c].replace(rbAtmDict.get(atomData[r][ci+8],''),'')
    12781379                if 'Atoms' in data['Drawing']:
    1279                     DrawAtomsReplaceByID(data['Drawing'],atomData[r],ID)
     1380                    DrawAtomsReplaceByID(data['Drawing'],ci+8,atomData[r],ID)
    12801381                wx.CallAfter(Paint)
    12811382
     
    13001401                ID = atomData[r][ci+8]
    13011402                if 'Atoms' in data['Drawing']:
    1302                     DrawAtomsReplaceByID(data['Drawing'],atomData[r],ID)
     1403                    DrawAtomsReplaceByID(data['Drawing'],ci+8,atomData[r],ID)
    13031404                SetupGeneral()
    13041405            else:
     
    16741775            FillAtomsGrid(Atoms)
    16751776            ID = atomData[indx[0]][ci+8]
    1676             DrawAtomsReplaceByID(data['Drawing'],atomData[indx[0]],ID)
     1777            DrawAtomsReplaceByID(data['Drawing'],ci+8,atomData[indx[0]],ID)
    16771778            G2plt.PlotStructure(G2frame,data)
    16781779        event.StopPropagation()
    16791780           
    1680     def DrawAtomsReplaceByID(drawingData,atom,ID):
     1781    def DrawAtomsReplaceByID(drawingData,loc,atom,ID):
    16811782        IDs = [ID,]
    16821783        atomData = drawingData['Atoms']
    1683         indx = G2mth.FindAtomIndexByIDs(atomData,IDs)
     1784        indx = G2mth.FindAtomIndexByIDs(atomData,loc,IDs)
    16841785        for ind in indx:
    16851786            atomData[ind] = MakeDrawAtom(atom,atomData[ind])
     
    18181919                            for r in indx:
    18191920                                ID = atomData[r][ci+8]
    1820                                 DrawAtomsReplaceByID(data['Drawing'],atomData[r],ID)
     1921                                DrawAtomsReplaceByID(data['Drawing'],ci+8,atomData[r],ID)
    18211922                    FillAtomsGrid(Atoms)
    18221923                dlg.Destroy()
     
    19042005                        XYZ = XYZ+cent+Cell
    19052006                        if Force:
    1906                             XYZ = G2spc.MoveToUnitCell(XYZ)
     2007                            XYZ,cell = G2spc.MoveToUnitCell(XYZ)
     2008                            Cell += cell
    19072009                        if New:
    19082010                            atom = copy.copy(atomData[ind])
     
    26672769                            atomData[r][c] = parms
    26682770                            drawAtoms.SetCellValue(r,c,parms)
    2669                         FindBondsDraw()
     2771                        FindBondsDraw(data)
    26702772                        G2plt.PlotStructure(G2frame,data)
    26712773                    dlg.Destroy()
     
    27092811                if drawAtoms.GetColLabelValue(c) in ['Style','Label']:
    27102812                    atomData[r][c] = drawAtoms.GetCellValue(r,c)
    2711                     FindBondsDraw()
     2813                    FindBondsDraw(data)
    27122814                elif drawAtoms.GetColLabelValue(c) == 'Color':
    27132815                    dlg = wx.ColourDialog(G2frame)
     
    27732875            labelChoice = [' ','type','name','number','residue','1-letter','chain']
    27742876            Types[9] = wg.GRID_VALUE_CHOICE+": ,type,name,number,residue,1-letter,chain"
    2775 #        elif generalData['Type'] == 'modulated':
    2776 #            Types += []
    2777 #            colLabels += []
    27782877        table = []
    27792878        rowLabels = []
     
    27822881            rowLabels.append(str(i))
    27832882
    2784         atomTable = G2G.Table(table,rowLabels=rowLabels,colLabels=colLabels,types=Types)
    2785         drawAtoms.SetTable(atomTable, True)
     2883        G2frame.atomTable = G2G.Table(table,rowLabels=rowLabels,colLabels=colLabels,types=Types)
     2884        drawAtoms.SetTable(G2frame.atomTable, True)
    27862885        drawAtoms.SetMargins(0,0)
    27872886        drawAtoms.AutoSizeColumns(True)
     
    28112910        G2frame.dataFrame.setSizePosLeft([600,300])
    28122911
    2813         FindBondsDraw()
     2912        FindBondsDraw(data)
    28142913        drawAtoms.ClearSelection()
    28152914#        G2plt.PlotStructure(G2frame,data)
     
    28332932                    drawAtoms.SetCellValue(r,cs,parms)
    28342933            dlg.Destroy()
    2835             FindBondsDraw()
     2934            FindBondsDraw(data)
    28362935            drawAtoms.ClearSelection()
    28372936            G2plt.PlotStructure(G2frame,data)
     
    29923091                        XYZ = XYZ+cent+Cell
    29933092                        if Force:
    2994                             XYZ = G2spc.MoveToUnitCell(XYZ)
     3093                            XYZ,cell = G2spc.MoveToUnitCell(XYZ)
     3094                            Cell += cell
    29953095                        atom = atomData[ind]
    29963096                        atom[cx:cx+3] = XYZ
     
    30853185                            atom[cuij:cuij+6] = item[1]
    30863186                            Opp = G2spc.Opposite(item[0])
    3087                             for xyz in Opp:
    3088                                 if noDuplicate(xyz,atomData):
    3089                                     cell = np.asarray(np.rint(xyz-atom[cx:cx+3]),dtype=np.int32)
    3090                                     cell = '1'+'+'+ \
    3091                                         str(cell[0])+','+str(cell[1])+','+str(cell[2])
    3092                                     atom[cx:cx+3] = xyz
    3093                                     atom[cx+3] = G2spc.StringOpsProd(cell,atom[cx+3],SGData)
     3187                            for key in Opp:
     3188                                if noDuplicate(Opp[key],atomData):
     3189                                    unit = np.array([eval(i) for i in key.split(',')])-item[3]
     3190                                    cell = '%d+%d,%d,%d'%(item[2],unit[0],unit[1],unit[2])
     3191                                    atom[cx:cx+3] = Opp[key]
     3192                                    atom[cx+3] = cell
    30943193                                    atomData.append(atom[:])
    30953194                    else:
     
    31013200                                +str(item[2][0])+','+str(item[2][1])+','+str(item[2][2])
    31023201                            Opp = G2spc.Opposite(item[0])
    3103                             for xyz in Opp:
    3104                                 if noDuplicate(xyz,atomData):
    3105                                     cell = np.asarray(np.rint(xyz-atom[cx:cx+3]),dtype=np.int32)
    3106                                     cell = '1'+'+'+ \
    3107                                         str(cell[0])+','+str(cell[1])+','+str(cell[2])
    3108                                     atom[cx:cx+3] = xyz
    3109                                     atom[cx+3] = G2spc.StringOpsProd(cell,atom[cx+3],SGData)
     3202                            for key in Opp:
     3203                                if noDuplicate(Opp[key],atomData):
     3204                                    unit = np.array([eval(i) for i in key.split(',')])-item[2]
     3205                                    cell = '%d+%d,%d,%d'%(item[1],unit[0],unit[1],unit[2])
     3206                                    atom[cx:cx+3] = Opp[key]
     3207                                    atom[cx+3] = cell
    31103208                                    atomData.append(atom[:])               
    31113209                    data['Drawing']['Atoms'] = atomData
     
    31163214            G2plt.PlotStructure(G2frame,data)
    31173215           
    3118     def FindBondsToo():                         #works but slow for large structures - keep as reference
    3119         cx,ct,cs,ci = data['Drawing']['atomPtrs']
    3120         atomData = data['Drawing']['Atoms']
    3121         generalData = data['General']
    3122         Amat,Bmat = G2lat.cell2AB(generalData['Cell'][1:7])
    3123         radii = generalData['BondRadii']
    3124         atomTypes = generalData['AtomTypes']
    3125         try:
    3126             indH = atomTypes.index('H')
    3127             radii[indH] = 0.5
    3128         except:
    3129             pass           
    3130         for atom in atomData:
    3131             atom[-1] = []
    3132         Atoms = []
    3133         for i,atom in enumerate(atomData):
    3134             Atoms.append([i,np.array(atom[cx:cx+3]),atom[cs],radii[atomTypes.index(atom[ct])]])
    3135         for atomA in Atoms:
    3136             if atomA[2] in ['lines','sticks','ellipsoids','balls & sticks','polyhedra']:
    3137                 for atomB in Atoms:                   
    3138                     Dx = atomB[1]-atomA[1]
    3139                     DX = np.inner(Amat,Dx)
    3140                     dist = np.sqrt(np.sum(DX**2))
    3141                     sumR = atomA[3]+atomB[3]
    3142                     if 0.5 < dist <= 0.85*sumR:
    3143                         i = atomA[0]
    3144                         if atomA[2] == 'polyhedra':
    3145                             atomData[i][-1].append(DX)
    3146                         elif atomB[1] != 'polyhedra':
    3147                             j = atomB[0]
    3148                             atomData[i][-1].append(Dx*atomA[3]/sumR)
    3149                             atomData[j][-1].append(-Dx*atomB[3]/sumR)
    3150                    
    3151     def FindBondsDraw():                    #uses numpy & masks - very fast even for proteins!
    3152         import numpy.ma as ma
    3153         cx,ct,cs,ci = data['Drawing']['atomPtrs']
    3154         hydro = data['Drawing']['showHydrogen']
    3155         atomData = data['Drawing']['Atoms']
    3156         generalData = data['General']
    3157         Amat,Bmat = G2lat.cell2AB(generalData['Cell'][1:7])
    3158         radii = generalData['BondRadii']
    3159         atomTypes = generalData['AtomTypes']
    3160         try:
    3161             indH = atomTypes.index('H')
    3162             radii[indH] = 0.5
    3163         except:
    3164             pass           
    3165         for atom in atomData:
    3166             atom[-2] = []               #clear out old bonds/polyhedra
    3167             atom[-1] = []
    3168         Indx = range(len(atomData))
    3169         Atoms = []
    3170         Styles = []
    3171         Radii = []
    3172         for atom in atomData:
    3173             Atoms.append(np.array(atom[cx:cx+3]))
    3174             Styles.append(atom[cs])
    3175             try:
    3176                 if not hydro and atom[ct] == 'H':
    3177                     Radii.append(0.0)
    3178                 else:
    3179                     Radii.append(radii[atomTypes.index(atom[ct])])
    3180             except ValueError:          #changed atom type!
    3181                 Radii.append(0.20)
    3182         Atoms = np.array(Atoms)
    3183         Radii = np.array(Radii)
    3184         IASR = zip(Indx,Atoms,Styles,Radii)
    3185         for atomA in IASR:
    3186             if atomA[2] in ['lines','sticks','ellipsoids','balls & sticks','polyhedra']:
    3187                 Dx = Atoms-atomA[1]
    3188                 dist = ma.masked_less(np.sqrt(np.sum(np.inner(Amat,Dx)**2,axis=0)),0.5) #gets rid of G2frame & disorder "bonds" < 0.5A
    3189                 sumR = atomA[3]+Radii
    3190                 IndB = ma.nonzero(ma.masked_greater(dist-data['Drawing']['radiusFactor']*sumR,0.))                 #get indices of bonded atoms
    3191                 i = atomA[0]
    3192                 for j in IndB[0]:
    3193                     if Styles[i] == 'polyhedra':
    3194                         atomData[i][-2].append(np.inner(Amat,Dx[j]))
    3195                     elif Styles[j] != 'polyhedra' and j > i:
    3196                         atomData[i][-2].append(Dx[j]*Radii[i]/sumR[j])
    3197                         atomData[j][-2].append(-Dx[j]*Radii[j]/sumR[j])
    3198                 if Styles[i] == 'polyhedra':
    3199                     Bonds = atomData[i][-2]
    3200                     Faces = []
    3201                     if len(Bonds) > 2:
    3202                         FaceGen = G2lat.uniqueCombinations(Bonds,3)     #N.B. this is a generator
    3203                         for face in FaceGen:
    3204                             vol = nl.det(face)
    3205                             if abs(vol) > 1. or len(Bonds) == 3:
    3206                                 if vol < 0.:
    3207                                     face = [face[0],face[2],face[1]]
    3208                                 face = np.array(face)
    3209                                 if not np.array([np.array(nl.det(face-bond))+0.0001 < 0 for bond in Bonds]).any():
    3210                                     norm = np.cross(face[1]-face[0],face[2]-face[0])
    3211                                     norm /= np.sqrt(np.sum(norm**2))
    3212                                     Faces.append([face,norm])
    3213                         atomData[i][-1] = Faces
    3214                        
    32153216    def DrawAtomsDelete(event):   
    32163217        indx = drawAtoms.GetSelectedRows()
     
    32273228       
    32283229    def OnReloadDrawAtoms(event):
    3229         data['Drawing']['Atoms'] = []
     3230        atomData = data['Atoms']
     3231        cx,ct,cs,ci = data['General']['AtomPtrs']
     3232        for atom in atomData:
     3233            ID = atom[ci+8]
     3234            DrawAtomsReplaceByID(data['Drawing'],ci+8,atom,ID)
    32303235        UpdateDrawAtoms()
    32313236        drawAtoms.ClearSelection()
     
    32353240    def DrawAtomsDeleteByIDs(IDs):
    32363241        atomData = data['Drawing']['Atoms']
    3237         indx = G2mth.FindAtomIndexByIDs(atomData,IDs)
     3242        loc = data['Drawing']['atomPtrs'][3]+8
     3243        indx = G2mth.FindAtomIndexByIDs(atomData,loc,IDs)
    32383244        indx.reverse()
    32393245        for ind in indx:
     
    32493255        elif colName == 'I/A':
    32503256            col = cs
    3251         indx = G2mth.FindAtomIndexByIDs(atomData,IDs)
     3257        indx = G2mth.FindAtomIndexByIDs(atomData,ci+8,IDs)
    32523258        for ind in indx:
    32533259            atomData[ind][col] = value
     
    33163322            atom = atomDData[i]
    33173323            xyz.append([i,]+atom[cn:cn+2]+atom[cx:cx+4]) #also gets Sym Op
    3318             id = G2mth.FindAtomIndexByIDs(atomData,[atom[cid],],False)[0]
     3324            id = G2mth.FindAtomIndexByIDs(atomData,cid,[atom[cid],],False)[0]
    33193325            Oxyz.append([id,]+atomData[id][cx+1:cx+4])
    33203326        DATData['Datoms'] = xyz
     
    35023508            def OnShowHyd(event):
    35033509                drawingData['showHydrogen'] = showHydrogen.GetValue()
    3504                 FindBondsDraw()
     3510                FindBondsDraw(data)
    35053511                G2plt.PlotStructure(G2frame,data)
    35063512               
    35073513            def OnShowRB(event):
    35083514                drawingData['showRigidBodies'] = showRB.GetValue()
    3509                 FindBondsDraw()
     3515                FindBondsDraw(data)
    35103516                G2plt.PlotStructure(G2frame,data)
    35113517               
     
    36143620                drawingData['radiusFactor'] = value
    36153621                radFactor.SetValue("%.2f"%(value))
    3616                 FindBondsDraw()
     3622                FindBondsDraw(data)
    36173623                G2plt.PlotStructure(G2frame,data)
    36183624           
Note: See TracChangeset for help on using the changeset viewer.