Changeset 4079


Ignore:
Timestamp:
Aug 5, 2019 3:33:15 PM (2 years ago)
Author:
vondreele
Message:

modify SHAPES for better output to GSAS-II
begin plotting capability for SHAPES output

Location:
trunk
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/G2shapes.py

    r4078 r4079  
    10831083    print (aString)
    10841084   
    1085     aList_summary = []
    1086     aList_summary.append(version_aString)
    1087     aList_summary.append(str(localtime))
     1085#    aList_summary = []
     1086#    aList_summary.append(version_aString)
     1087#    aList_summary.append(str(localtime))
    10881088   
    10891089    ######################
     
    13771377                         
    13781378        aString  = 'Trial:' + file_no
    1379         aList_summary.append(aString)
     1379#        aList_summary.append(aString)
    13801380   
    13811381        file_beads = prefix + 'beads_' + file_no + '.pdb'
    1382         file_pr = prefix + 'pr_calc_' + file_no + '.dat'
     1382#        file_pr = prefix + 'pr_calc_' + file_no + '.dat'
    13831383        file_psv = prefix + 'psv_shape_' + file_no + '.pdb'
    13841384#        file_intensity = prefix + 'intensity_' + file_no + '.dat'
     
    17801780        print (aString)
    17811781   
     1782        Phases.append([file_beads.split('.')[0],aList_beads_x,aList_beads_y,aList_beads_z])
    17821783        # Write out beads as pseudo a PDB file
    1783         Phases.append({'pname':file_beads.split('.')[0],
    1784             'Atoms':[aList_beads_x,aList_beads_y,aList_beads_z,file_beads]})
    17851784        if pdbOut:
    17861785            pdb_writer(aList_beads_x,aList_beads_y,aList_beads_z,file_beads,aString)
     
    18301829        # Write input and model P(r)
    18311830#        pr_writer(aList_pr,aList_r,aList_pr_model,file_pr)
    1832         PRcalc.append([aList_r,aList_pr,aList_pr_model,])
     1831        PRcalc.append([aList_r,aList_pr,aList_pr_model,delta_hist_sum])
    18331832   
    18341833        # Calculate comparison versus intensities
     
    18481847   
    18491848            # Write output intensity file
    1850             Patterns.append([aList_q,aList_i,aList_i_calc,aString])
     1849            Patterns.append([aList_q,aList_i,aList_i_calc,rvalue])
    18511850#            write_all_data(file_intensity,aList_q,aList_i,aList_i_calc,aString)
    18521851   
     
    18691868                  ' PSV:' + str(fraction_psv)
    18701869   
    1871         Phases.append({'pname':file_psv.split('.')[0],
    1872             'Atoms':[aList_box_x,aList_box_y,aList_box_z]})
     1870        Phases.append([file_psv.split('.')[0],aList_box_x,aList_box_y,aList_box_z])
    18731871        if pdbOut:
    18741872            pdb_writer(aList_box_x,aList_box_y,aList_box_z,file_psv,aString)
     
    18781876        aString = 'P(r) dif:' + str(delta_hist_sum) + ' E:' \
    18791877                         + str(vdw_all) + ' CHISQ:' + str(chi_sq) + ' PSV:' + str(fraction_psv)
    1880         aList_summary.append(aString)               
     1878#        aList_summary.append(aString)               
    18811879   
    18821880        i_soln = i_soln + 1
     
    19011899    print (aString)
    19021900                         
    1903     aList_summary.append(str(localtime))
    1904    
    1905     # Create summary
    1906    
    1907     aFile_log = prefix + 'shapes_summary.log'
    1908     num_lines = len(aList_summary)
    1909    
    1910     file = open(aFile_log,'w')
    1911     i = 0
    1912     while i < num_lines:
    1913         aString = str(aList_summary[i])
    1914         file.write(aString)
    1915         file.write('\n')
    1916         i = i + 1
    1917     file.close()
     1901#    aList_summary.append(str(localtime))
     1902#   
     1903#    # Create summary
     1904#   
     1905#    aFile_log = prefix + 'shapes_summary.log'
     1906#    num_lines = len(aList_summary)
     1907#   
     1908#    file = open(aFile_log,'w')
     1909#    i = 0
     1910#    while i < num_lines:
     1911#        aString = str(aList_summary[i])
     1912#        file.write(aString)
     1913#        file.write('\n')
     1914#        i = i + 1
     1915#    file.close()
    19181916
    19191917    return Phases,Patterns,PRcalc
  • trunk/GSASIIplot.py

    r4077 r4079  
    46794679        Ymax = np.amax(XYlist.T[1])
    46804680        dy = 0.02*(Ymax-Ymin)
    4681         Plot.set_xlim(Xmin-dx,Xmax+dx)
    4682         Plot.set_ylim(Ymin-dy,Ymax+dy)
     4681        try:
     4682            Plot.set_xlim(Xmin-dx,Xmax+dx)
     4683            Plot.set_ylim(Ymin-dy,Ymax+dy)
     4684        except:
     4685            pass
    46834686        if Peaks == None:
    46844687            normcl = mpcls.Normalize(Ymin,Ymax)
     
    52255228        elif 'Pair' in PlotText:
    52265229            PlotSASDPairDist(G2frame)
     5230           
    52275231   
    52285232    def OnMotion(event):
     
    91349138
    91359139################################################################################
     9140#### Plot Bead Model
     9141################################################################################
     9142
     9143def PlotBeadModel(G2frame,Atoms,defaults):
     9144    '''Bead modelplotting package. For bead models from SHAPES
     9145    '''
     9146
     9147    Mydir = G2frame.dirname
     9148    Rd = np.array([255,0,0])
     9149    Gr = np.array([0,255,0])
     9150    Bl = np.array([0,0,255])
     9151    uBox = np.array([[0,0,0],[1,0,0],[0,1,0],[0,0,1]])
     9152    uEdges = np.array([[uBox[0],uBox[1]],[uBox[0],uBox[2]],[uBox[0],uBox[3]]])
     9153    uColors = [Rd,Gr,Bl]
     9154    XYZ = np.array(Atoms[1:]).T      #don't mess with original!
     9155
     9156#    def SetRBOrigin():
     9157#        page = getSelection()
     9158#        if page:
     9159#            if G2frame.GetPageText(page) == 'Rigid bodies':
     9160#                G2frame.MapPeaksTable.SetData(mapPeaks)
     9161#                panel = G2frame.GetPage(page).GetChildren()
     9162#                names = [child.GetName() for child in panel]
     9163#                panel[names.index('grid window')].Refresh()
     9164           
     9165    def OnMouseDown(event):
     9166        xy = event.GetPosition()
     9167        defaults['oldxy'] = list(xy)
     9168
     9169    def OnMouseMove(event):
     9170        newxy = event.GetPosition()
     9171                               
     9172        if event.Dragging():
     9173            if event.LeftIsDown():
     9174                SetRotation(newxy)
     9175                Q = defaults['Quaternion']
     9176                G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
     9177#            elif event.RightIsDown():
     9178#                SetRBOrigin(newxy)
     9179            elif event.MiddleIsDown():
     9180                SetRotationZ(newxy)
     9181                Q = defaults['Quaternion']
     9182                G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
     9183            Draw('move')
     9184       
     9185    def OnMouseWheel(event):
     9186        defaults['cameraPos'] += event.GetWheelRotation()/24
     9187        defaults['cameraPos'] = max(10,min(500,defaults['cameraPos']))
     9188        G2frame.G2plotNB.status.SetStatusText('New camera distance: %.2f'%(defaults['cameraPos']),1)
     9189        Draw('wheel')
     9190       
     9191    def SetBackground():
     9192        R,G,B,A = Page.camera['backColor']
     9193        GL.glClearColor(R,G,B,A)
     9194        GL.glClear(GL.GL_COLOR_BUFFER_BIT | GL.GL_DEPTH_BUFFER_BIT)
     9195       
     9196    def SetLights():
     9197        GL.glEnable(GL.GL_DEPTH_TEST)
     9198        GL.glShadeModel(GL.GL_FLAT)
     9199        GL.glEnable(GL.GL_LIGHTING)
     9200        GL.glEnable(GL.GL_LIGHT0)
     9201        GL.glLightModeli(GL.GL_LIGHT_MODEL_TWO_SIDE,0)
     9202        GL.glLightfv(GL.GL_LIGHT0,GL.GL_AMBIENT,[1,1,1,.8])
     9203        GL.glLightfv(GL.GL_LIGHT0,GL.GL_DIFFUSE,[1,1,1,1])
     9204       
     9205    def SetRotation(newxy):
     9206#first get rotation vector in screen coords. & angle increment       
     9207        oldxy = defaults['oldxy']
     9208        if not len(oldxy): oldxy = list(newxy)
     9209        dxy = newxy-oldxy
     9210        defaults['oldxy'] = list(newxy)
     9211        V = np.array([dxy[1],dxy[0],0.])
     9212        A = 0.25*np.sqrt(dxy[0]**2+dxy[1]**2)
     9213# next transform vector back to xtal coordinates via inverse quaternion
     9214# & make new quaternion
     9215        Q = defaults['Quaternion']
     9216        V = G2mth.prodQVQ(G2mth.invQ(Q),V)
     9217        DQ = G2mth.AVdeg2Q(A,V)
     9218        Q = G2mth.prodQQ(Q,DQ)
     9219        defaults['Quaternion'] = Q
     9220# finally get new view vector - last row of rotation matrix
     9221        VD = G2mth.Q2Mat(Q)[2]
     9222        VD /= np.sqrt(np.sum(VD**2))
     9223        defaults['viewDir'] = VD
     9224       
     9225    def SetRotationZ(newxy):                       
     9226#first get rotation vector (= view vector) in screen coords. & angle increment       
     9227        View = GL.glGetIntegerv(GL.GL_VIEWPORT)
     9228        cent = [View[2]/2,View[3]/2]
     9229        oldxy = defaults['oldxy']
     9230        if not len(oldxy): oldxy = list(newxy)
     9231        dxy = newxy-oldxy
     9232        defaults['oldxy'] = list(newxy)
     9233        V = defaults['viewDir']
     9234        A = [0,0]
     9235        A[0] = dxy[1]*.25
     9236        A[1] = dxy[0]*.25
     9237        if newxy[0] > cent[0]:
     9238            A[0] *= -1
     9239        if newxy[1] < cent[1]:
     9240            A[1] *= -1       
     9241# next transform vector back to xtal coordinates & make new quaternion
     9242        Q = defaults['Quaternion']
     9243        Qx = G2mth.AVdeg2Q(A[0],V)
     9244        Qy = G2mth.AVdeg2Q(A[1],V)
     9245        Q = G2mth.prodQQ(Q,Qx)
     9246        Q = G2mth.prodQQ(Q,Qy)
     9247        defaults['Quaternion'] = Q
     9248
     9249    def RenderUnitVectors(x,y,z):
     9250        GL.glEnable(GL.GL_COLOR_MATERIAL)
     9251        GL.glLineWidth(1)
     9252        GL.glPushMatrix()
     9253        GL.glTranslate(x,y,z)
     9254        GL.glBegin(GL.GL_LINES)
     9255        for line,color in zip(uEdges,uColors):
     9256            GL.glColor3ubv(color)
     9257            GL.glVertex3fv(-line[1])
     9258            GL.glVertex3fv(line[1])
     9259        GL.glEnd()
     9260        GL.glPopMatrix()
     9261        GL.glColor4ubv([0,0,0,0])
     9262        GL.glDisable(GL.GL_COLOR_MATERIAL)
     9263               
     9264    def RenderSphere(x,y,z,radius,color):
     9265        GL.glMaterialfv(GL.GL_FRONT_AND_BACK,GL.GL_DIFFUSE,color)
     9266        GL.glPushMatrix()
     9267        GL.glTranslate(x,y,z)
     9268        q = GLU.gluNewQuadric()
     9269        GLU.gluSphere(q,radius,20,10)
     9270        GL.glPopMatrix()
     9271       
     9272    def Draw(caller=''):
     9273#useful debug?       
     9274#        if caller:
     9275#            print caller
     9276# end of useful debug
     9277        cPos = defaults['cameraPos']
     9278        VS = np.array(Page.canvas.GetSize())
     9279        aspect = float(VS[0])/float(VS[1])
     9280        Q = defaults['Quaternion']
     9281        SetBackground()
     9282        GL.glInitNames()
     9283        GL.glPushName(0)
     9284       
     9285        GL.glMatrixMode(GL.GL_PROJECTION)
     9286        GL.glLoadIdentity()
     9287        GL.glViewport(0,0,VS[0],VS[1])
     9288        GLU.gluPerspective(20.,aspect,1.,500.)
     9289        GLU.gluLookAt(0,0,cPos,0,0,0,0,1,0)
     9290        SetLights()           
     9291           
     9292        GL.glMatrixMode(GL.GL_MODELVIEW)
     9293        GL.glLoadIdentity()
     9294        matRot = G2mth.Q2Mat(Q)
     9295        matRot = np.concatenate((np.concatenate((matRot,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
     9296        GL.glMultMatrixf(matRot.T)
     9297        RenderUnitVectors(0.,0.,0.)
     9298        radius = 2.0
     9299        for iat,atom in enumerate(XYZ):
     9300            x,y,z = atom
     9301            color = (144,144,144)
     9302            RenderSphere(x,y,z,radius,color)
     9303        try:
     9304            if Page.context: Page.canvas.SetCurrent(Page.context)
     9305        except:
     9306            pass
     9307        Page.canvas.SwapBuffers()
     9308
     9309    def OnSize(event):
     9310        Draw('size')
     9311       
     9312    def OnFocus(event):
     9313        Draw('focus')
     9314       
     9315    def OnKeyBox(event):
     9316        mode = cb.GetValue()
     9317        if mode in ['jpeg','bmp','tiff',]:
     9318            try:
     9319                import Image as Im
     9320            except ImportError:
     9321                try:
     9322                    from PIL import Image as Im
     9323                except ImportError:
     9324                    print ("PIL/pillow Image module not present. Cannot save images without this")
     9325                    raise Exception("PIL/pillow Image module not found")
     9326           
     9327            Fname = os.path.join(Mydir,Page.name+'.'+mode)
     9328            print (Fname+' saved')
     9329            size = Page.canvas.GetSize()
     9330            GL.glPixelStorei(GL.GL_UNPACK_ALIGNMENT, 1)
     9331            Pix = GL.glReadPixels(0,0,size[0],size[1],GL.GL_RGB, GL.GL_UNSIGNED_BYTE)
     9332            im = Im.new("RGB", (size[0],size[1]))
     9333            try:
     9334                im.frombytes(Pix)
     9335            except AttributeError:
     9336                im.fromstring(Pix)
     9337            im = im.transpose(Im.FLIP_TOP_BOTTOM)
     9338            im.save(Fname,mode)
     9339            cb.SetValue(' save as/key:')
     9340            G2frame.G2plotNB.status.SetStatusText('Drawing saved to: '+Fname,1)
     9341
     9342    # PlotRigidBody execution starts here (N.B. initialization above)
     9343    new,plotNum,Page,Plot,lim = G2frame.G2plotNB.FindPlotTab('Bead model','ogl')
     9344    if new:
     9345        Page.views = False
     9346    Page.name = Atoms[0]
     9347    Page.Choice = None
     9348    choice = [' save as:','jpeg','tiff','bmp',]
     9349    cb = wx.ComboBox(G2frame.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,choices=choice)
     9350    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
     9351    cb.SetValue(' save as/key:')
     9352   
     9353    Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
     9354    Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
     9355    Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
     9356    Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
     9357    Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
     9358    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
     9359    Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
     9360    Page.camera['position'] = defaults['cameraPos']
     9361    Page.camera['backColor'] = np.array([0,0,0,0])
     9362    try:
     9363        if Page.context: Page.canvas.SetCurrent(Page.context)
     9364    except:
     9365        pass
     9366    Draw('main')
     9367    Draw('main')    #to fill both buffers so save works
     9368
     9369################################################################################
    91369370#### Plot Rigid Body
    91379371################################################################################
     
    93779611        Draw('size')
    93789612       
     9613    def OnFocus(event):
     9614        Draw('focus')
     9615       
    93799616    def OnKeyBox(event):
    93809617        mode = cb.GetValue()
     
    94229659    Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
    94239660    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
     9661    Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
    94249662    Page.camera['position'] = defaults['cameraPos']
    94259663    Page.camera['backColor'] = np.array([0,0,0,0])
     
    983110069        Draw('size')
    983210070       
     10071    def OnFocus(event):
     10072        Draw('focus')
     10073       
    983310074    # PlotLayers execution starts here
    983410075    cell = Layers['Cell'][1:7]
     
    987610117    Page.canvas.Bind(wx.EVT_KEY_UP, OnPlotKeyPress)
    987710118    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
     10119    Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
    987810120    Page.camera['position'] = defaults['cameraPos']
    987910121    Page.camera['backColor'] = np.array([0,0,0,0])
  • trunk/GSASIIpwdGUI.py

    r4078 r4079  
    358358        'Pair':{'Method':'Regularization','MaxRadius':100.,'NBins':100,'Errors':'User',
    359359                'Percent error':2.5,'Background':[0,False],'Distribution':[],
    360                 'Moore':20,'Dist G':100.,},           
     360                'Moore':20,'Dist G':100.,'Result':[],},           
    361361        'Particle':{'Matrix':{'Name':'vacuum','VolFrac':[0.0,False]},'Levels':[],},
    362362        'Shapes':{'outName':'run','NumAA':100,'Niter':1,'AAscale':1.0,'Symm':1,'bias-z':0.0,
     
    52245224        data['BackFile'] = ''
    52255225    if 'Pair' not in data:
    5226         data['Pair'] = {'Method':'Regularization','MaxRadius':100.,'NBins':100,'Errors':'User',
     5226        data['Pair'] = {'Method':'Regularization','MaxRadius':100.,'NBins':100,'Errors':'User','Result':[],
    52275227            'Percent error':2.5,'Background':[0,False],'Distribution':[],'Moore':[20,False],'Dist G':100.,} 
    52285228    if 'Shapes' not in data:
    52295229        data['Shapes'] = {'outName':'run','NumAA':100,'Niter':1,'AAscale':1.0,'Symm':1,'bias-z':0.0,
    52305230            'inflateV':1.0,'AAglue':0.0,'pdbOut':False}
     5231    plotDefaults = {'oldxy':[0.,0.],'Quaternion':[0.,0.,0.,1.],'cameraPos':150.,'viewDir':[0,0,1],}
    52315232
    52325233    #end patches
     
    54285429      doi: xxx''',
    54295430      caption='Program Shapes',style=wx.ICON_INFORMATION)
    5430             Phases,Patterns,PRcalc = G2shapes.G2shapes(Profile,ProfDict,Limits,data)
    5431            
    5432            
     5431            data['Pair']['Result'] = []       #clear old results (if any) for now
     5432            data['Pair']['Result'] = G2shapes.G2shapes(Profile,ProfDict,Limits,data)
     5433            wx.CallAfter(UpdateModelsGrid,G2frame,data)
    54335434           
    54345435    def OnUnDo(event):
     
    56545655        def OnPDBout(event):
    56555656            data['Shapes']['pdbOut'] = not data['Shapes']['pdbOut']
     5657           
     5658        def OnShapeSelect(event):
     5659            r,c =  event.GetRow(),event.GetCol()
     5660            shapeTable.SetValue(r,c,False)
     5661            selAtoms = Atoms[2*r+(c-1)]
     5662            pattern = Patterns[r]
     5663            PRvals = PRcalc[r]
     5664            print('%s %d'%('num. beads',len(selAtoms[1])))
     5665            print('%s %.3f'%('selected r value',pattern[-1]))
     5666            print('%s %.3f'%('selected Delta P(r)',PRvals[-1]))
     5667            G2plt.PlotSASDPairDist(G2frame)
     5668            RefreshPlots(True)
     5669           
     5670            G2plt.PlotBeadModel(G2frame,selAtoms,plotDefaults)
     5671           
     5672           
    56565673       
    56575674        shapeSizer = wx.BoxSizer(wx.VERTICAL)
     
    56935710       
    56945711        shapeSizer.Add(parmSizer)
     5712       
     5713        if len(data['Pair'].get('Result',[])):
     5714            shapeSizer.Add(wx.StaticText(G2frame.dataWindow,label=' SHAPES run results:'),0,WACV)
     5715            Atoms,Patterns,PRcalc = data['Pair']['Result']
     5716            colLabels = ['name','show beads','show shape','Rvalue','P(r) dif','Nbeads','Nshape']
     5717            Types = [wg.GRID_VALUE_STRING,]+2*[wg.GRID_VALUE_BOOL,]+2*[wg.GRID_VALUE_FLOAT+':10,3',]+2*[wg.GRID_VALUE_LONG,]
     5718            rowLabels = [str(i) for i in range(len(Patterns))]
     5719            tableVals = []
     5720            for i in range(len(Patterns)):
     5721                tableVals.append([Atoms[2*i][0],False,False,Patterns[i][-1],PRcalc[i][-1],len(Atoms[2*1][1]),len(Atoms[2*i+1][1])])
     5722            shapeTable = G2G.Table(tableVals,rowLabels=rowLabels,colLabels=colLabels,types=Types)
     5723            ShapesResult = G2G.GSGrid(G2frame.dataWindow)
     5724            ShapesResult.SetTable(shapeTable,True)
     5725            ShapesResult.AutoSizeColumns(False)
     5726            ShapesResult.Bind(wg.EVT_GRID_CELL_CHANGED, OnShapeSelect)
     5727            for r in range(len(Patterns)):
     5728                for c in range(7):
     5729                    if c in [1,2]:
     5730                        ShapesResult.SetReadOnly(r,c,isReadOnly=False)
     5731                    else:
     5732                        ShapesResult.SetReadOnly(r,c,isReadOnly=True)
     5733   
     5734            shapeSizer.Add(ShapesResult,0,WACV)
    56955735        return shapeSizer
    56965736
Note: See TracChangeset for help on using the changeset viewer.