Changeset 5520

Timestamp:

Mar 23, 2023 1:03:34 PM (2 years ago)

Author:

vondreele

Message:

New version of spinning RBs - radius is now a parameter of each shell. all appear to work

Location:

trunk

Files:

• GSASIIconstrGUI.py (modified) (5 diffs)
• GSASIIlattice.py (modified) (1 diff)
• GSASIIobj.py (modified) (1 diff)
• GSASIIphsGUI.py (modified) (8 diffs)
• GSASIIplot.py (modified) (1 diff)
• GSASIIstrIO.py (modified) (10 diffs)
• GSASIIstrMain.py (modified) (1 diff)
• GSASIIstrMath.py (modified) (8 diffs)

trunk/GSASIIconstrGUI.py ¶

@@ -3027 +3027 @@
         
         rbid = ran.randint(0,sys.maxsize)
-        radius = 1.0
         atType = 'C'
         rbType = 'Q'
@@ -3035 +3034 @@
         name = G2obj.MakeUniqueLabel(name,namelist)
         atColor = G2elem.GetAtomInfo(atType)['Color']
-        data['Spin'][rbid] = {'RBname':name,'Natoms':Natoms,'radius':[radius,False],'atType':atType,'rbType':rbType,
-            'useCount':0,'nSH':0,'SHC':[{},],'Matrix':np.eye(3),'rbPos':np.zeros(3),'atColor':atColor}
+        data['Spin'][rbid] = {'RBname':name,'Natoms':Natoms,'atType':atType,'rbType':rbType,
+            'useCount':0,'nSH':0,'SHC':[{},],'Radius':[1.0,False],'Matrix':np.eye(3),'rbPos':np.zeros(3),'atColor':atColor}
         data['RBIds']['Spin'].append(rbid)
         UpdateSpinRB()
@@ -3404 +3403 @@
                     Obj.ChangeValue(El)
                    
-        def OnRefSel(event):
-            ObjId = event.GetId()
-            data['Spin'][Indx[ObjId]]['radius'][1] = not data['Spin'][Indx[ObjId]]['radius'][1]
-
         def OnSymSel(event):
             ObjId = event.GetId()
@@ -3434 +3429 @@
         Indx = {}
         SpinRBSizer.Add(wx.StaticText(SpinRBDisplay,label=' Spinning rigid body shells/nonspherical atoms (radius=0):'),0,WACV)
-        bodSizer = wx.FlexGridSizer(0,7,5,5)
-        for item in ['Name','Type','RB sym','Atom','Number','radius','refine']:
+        bodSizer = wx.FlexGridSizer(0,5,5,5)
+        for item in ['Name','Type','RB sym','Atom','Number']:
             bodSizer.Add(wx.StaticText(SpinRBDisplay,label=item))
         for ibod,spinID in enumerate(data['Spin']):
@@ -3453 +3448 @@
             bodSizer.Add(atSel,0)
             bodSizer.Add(G2G.ValidatedTxtCtrl(SpinRBDisplay,data['Spin'][spinID],'Natoms'))
-            bodSizer.Add(G2G.ValidatedTxtCtrl(SpinRBDisplay,data['Spin'][spinID]['radius'],0))
-            refSel = wx.CheckBox(SpinRBDisplay)
-            refSel.SetValue(data['Spin'][spinID]['radius'][1])
-            Indx[refSel.GetId()] = spinID
-            refSel.Bind(wx.EVT_CHECKBOX,OnRefSel)
-            bodSizer.Add(refSel,0)
         
         SpinRBSizer.Add(bodSizer)

trunk/GSASIIlattice.py ¶

@@ -2706 +2706 @@
 def SHarmcal(SytSym,SHFln,psi,gam):
     '''Perform a surface spherical harmonics computation.
+    Presently only used for plotting
     Note that the the number of gam values must either be 1 or must match psi
     

trunk/GSASIIobj.py ¶

@@ -2041 +2041 @@
         'RBSO([aijk])' : 'Spinning rigid body orientation parameter \\1',
         'RBSP([xyz])' : 'Spinning rigid body \\1 position parameter',
+        'RBSShRadius' : 'Spinning rigid body radius',
         'RBSShC([1-20,1-20])'  : 'Spinning rigid body sph. harmonics term',
         'constr([0-9]*)' : 'Generated degree of freedom from constraint',

trunk/GSASIIphsGUI.py ¶

@@ -11581 +11581 @@
                 RBObj['RBId'].append(rbIds[selection])
                 RBObj['SHC'].append({})
-                for name in ['atColor','atType','Natoms','nSH','radius','RBname','RBsym']:
+                for name in ['atColor','atType','Natoms','nSH','RBname','RBsym','Radius']:
                     RBObj[name].append(rbData[name])
                 RBObj['hide'].append(False)
@@ -11612 +11612 @@
                     iSh,name = Indx[Obj.GetId()]
                     RBObj['SHC'][iSh][name][2] = not RBObj['SHC'][iSh][name][2]
+                    
+                def OnRadRef(event):
+                    Obj = event.GetEventObject()
+                    iSh = Indx[Obj.GetId()]
+                    RBObj['Radius'][iSh][1] = not RBObj['Radius'][iSh][1]
                 
                 def NewSHC(invalid,value,tc):
                     G2plt.PlotStructure(G2frame,data)
-                    
+                                       
                 def OnDelShell(event):
                     Obj = event.GetEventObject()
@@ -11622 +11627 @@
                     RBData['Spin'][rbId]['useCount'] -= 1
                     RBData['Spin'][rbId]['useCount'] = max(0,RBData['Spin'][rbId]['useCount'])
-                    for name in ['atColor','atType','Natoms','nSH','radius','RBId','RBname','RBsym','SHC']:
+                    for name in ['atColor','atType','Natoms','nSH','RBId','RBname','RBsym','SHC']:
                         del RBObj[name][iSh]
                     G2plt.PlotStructure(G2frame,data)
@@ -11630 +11635 @@
                 shSizer = wx.BoxSizer(wx.VERTICAL)
                 for iSh,nSh in enumerate(RBObj['nSH']):
+                    #patch
+                    if 'Radius' not in RBObj:
+                        RBObj['Radius'] = [[1.0,False] for i in range(len(RBObj['nSH']))]
+                    #end patch
                     rbId = RBObj['RBId'][iSh]
                     RBObj['atType'][iSh] = RBData['Spin'][rbId]['atType']
@@ -11654 +11663 @@
                     shoSizer.Add(shOrder,0,WACV)
                     if RBObj['nSH'][iSh]>0:
-                        shoSizer.Add(wx.StaticText(RigidBodies,label=" 'c' for cubic harmonic term"),0,WACV)
+                        shoSizer.Add(wx.StaticText(RigidBodies,label=" 'c' for cubic harmonic term. "),0,WACV)
+                    shoSizer.Add(wx.StaticText(RigidBodies,label=' Radius: '),0,WACV)
+                    shoSizer.Add(G2G.ValidatedTxtCtrl(RigidBodies,RBObj['Radius'][iSh],0,nDig=(8,5),xmin=0.0,xmax=5.0,
+                        typeHint=float,size=(70,-1),OnLeave=NewSHC),0,WACV)
+                    radref = wx.CheckBox(RigidBodies,label=' refine? ')
+                    radref.SetValue(RBObj['Radius'][iSh][1])
+                    radref.Bind(wx.EVT_CHECKBOX,OnRadRef)
+                    Indx[radref.GetId()] = iSh
+                    shoSizer.Add(radref,0,WACV)
                     shSizer.Add(shoSizer)
                     if not RBObj['nSH'][iSh]:
@@ -12287 +12304 @@
                     data['RBModels'][rbType] = []
                 if rbType == 'Spin':    #convert items to lists of shells
-                    for name in ['atColor','atType','Natoms','nSH','radius','RBId','RBname','RBsym']:
+                    for name in ['atColor','atType','Natoms','nSH','Radius','RBId','RBname','RBsym']:
                         item = rbObj[name]                        
                         rbObj[name] = [item,] 
@@ -12942 +12959 @@
         if rbType == 'Spin':
             data['testRBObj']['rbAtTypes'] = [RBData[rbType][rbId]['rbType'],] 
-            data['testRBObj']['AtInfo'] = {RBData[rbType][rbId]['rbType']:[RBData[rbType][rbId]['radius'],(128, 128, 255)],}
+            data['testRBObj']['AtInfo'] = {RBData[rbType][rbId]['rbType']:[RBData[rbType][rbId]['Radius'],(128, 128, 255)],}
             data['testRBObj']['rbType'] = rbType
             data['testRBObj']['rbData'] = RBData
@@ -12975 +12992 @@
                     'Orient':[[0.,0.,0.,1.],' '],'Ids':[],'RBId':rbId,'Torsions':[],
                     'numChain':'','RBname':RBData[rbType][rbId]['RBname']}
+        # if rbType == 'Spin':
+        #     data['testRBObj']['rbObj']['Radius'] = [1.0,False]
         data['testRBObj']['torAtms'] = []                
         for item in RBData[rbType][rbId].get('rbSeq',[]):

trunk/GSASIIplot.py ¶

@@ -10349 +10349 @@
                     if SpnData is not None:
                         SytSym = G2spc.SytSym(atom[cx:cx+3],SGData)[0]
-                        radius = SpnData['radius']
+                        radius = SpnData['Radius']
                         atColor = SpnData['atColor']
                         symAxis = np.array(SpnData.get('symAxis',[0,0,1]))

trunk/GSASIIstrIO.py ¶

@@ -1059 +1059 @@
         pFile.write('Spinning RB name: %s atom type: %s, no.atoms: %d, No. times used: %d\n'%
             (RBModel['RBname'],RBModel['atType'],RBModel['Natoms'],RBModel['useCount']))
-        for i in WriteSpnRBModel(RBModel):
-            pFile.write(i)
 
     def PrintResRBModel(RBModel):
@@ -1083 +1081 @@
         for irb,item in enumerate(rbIds['Spin']):
             if rigidbodyDict['Spin'][item]['useCount']:
-                RBradius = rigidbodyDict['Spin'][item]['radius']
-                pid = '::RBS;0:'+str(irb)
-                rbDict[pid] = RBradius[0]
-                if RBradius[1]:
-                    rbVary.append(pid)
                 if Print:
                     pFile.write('\nSpinning rigid body model:\n')
@@ -1132 +1125 @@
         pFile.write(sigstr+'\n')
 
-    def PrintRBSpnandSig(SpinRB,SpinSig):
-        pFile.write('\n Spinning radius for %s:\n'%SpinRB['RBname'])
-        namstr = '  names :'
-        valstr = '  values:'
-        sigstr = '  esds  :'
-        for i,[val,sig] in enumerate(zip(SpinRB['radius'],SpinSig)):
-            namstr += '%12s'%('Spin '+str(i))
-            valstr += '%12.4f'%(val)
-            if sig:
-                sigstr += '%12.4f'%(sig)
-            else:
-                sigstr += 12*' '
-        pFile.write(namstr+'\n')
-        pFile.write(valstr+'\n')
-        pFile.write(sigstr+'\n')        
-        
     RBIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[],'Spin':[]})  #these are lists of rbIds
-    if not RBIds['Vector'] and not RBIds['Spin']:
+    if not RBIds['Vector']:
         return
     for irb,item in enumerate(RBIds['Vector']):
@@ -1160 +1137 @@
                     VectSig.append(sigDict[name])
             PrintRBVectandSig(rigidbodyDict['Vector'][item],VectSig)    
-    for irb,item in enumerate(RBIds['Spin']):
-        if rigidbodyDict['Spin'][item]['useCount']:
-            SpinSig = []
-            name = '::RBS;'+str(0)+':'+str(irb)
-            if name in sigDict:
-                SpinSig.append(sigDict[name])
-                
-            
-            PrintRBSpnandSig(rigidbodyDict['Spin'][item],SpinSig)    
         
 ################################################################################
@@ -1343 +1311 @@
                         iFin = min(iBeg+6,nCoeff)
                         for block in range(nBlock):
-                            ptlbls = ' names :'
-                            ptstr =  ' values:'
-                            ptref =  ' refine:'
+                            if not block:
+                                ptlbls = ' names :%12s'%'Radius'
+                                ptstr =  ' values:%12.4f'%RB['Radius'][ish][0]
+                                ptref =  ' refine:%12s'%RB['Radius'][ish][1]                               
+                            else:
+                                ptlbls = ' names :'
+                                ptstr =  ' values:'
+                                ptref =  ' refine:'
                             for item in SHkeys[iBeg:iFin]:
                                 ptlbls += '%12s'%(item)
@@ -1592 +1565 @@
                 continue
             rbid = str(rbids.index(RB['RBId'][ish]))
+            name = '%sRBSSh;%d;Radius:%s:%s'%(pfx,ish,iAt,rbid)
+            phaseDict[name] = RB['Radius'][ish][0]
+            if RB['Radius'][ish][1]:
+                phaseVary += [name,]
             pfxRB = '%sRBSSh;%d;'%(pfx,ish)
             for i,shcof in enumerate(Shcof):
@@ -2569 +2546 @@
                         jrb = SRBIds.index(RBObj['RBId'][ish])  #spin rb no.
                         rbsx = ':%d:%d'%(iAt,jrb)
-                        pFile.write(' Spin rigid body parameters:\n')
-                        PrintRBObjPOAndSig('RBS',rbsx)
+                        pFile.write(' Spin rigid body parameters for at. no. %d; shell %d:\n'%(iAt,ish))
+                        if not ish:
+                            PrintRBObjPOAndSig('RBS',rbsx[1:])  #skip leading ':'
                         PrintRBObjSHCAndSig('RBSSh;%d;'%ish,RBObj['SHC'][ish],rbsx)
             atomsSig = {}
@@ -4202 +4180 @@
         if name in sigDict:
             sigstr += '%12.4f'%(sigDict[name])
+        else:
+            sigstr += 12*' '
     out.append(namstr+'\n')
     out.append(valstr+'\n')
@@ -4212 +4192 @@
     '''
     out = []
-    namstr = '  names :'
-    valstr = '  values:'
+    name = pfx+rbfx+'Radius'+rbsx
+    namstr = '  names :%12s'%'Radius'
+    valstr = '  values:%12.4f'%parmDict[name]
     sigstr = '  esds  :'
+    if name in sigDict:
+        sigstr += '%12.4f'%sigDict[name]
+        [name]
+    else:
+        sigstr += 12*' '
     out.append(' Sp.harm.:\n')
     for item in SHC:
@@ -4222 +4208 @@
         if name in sigDict:
             sigstr += '%12.4f'%(sigDict[name])
+        else:
+            sigstr += 12*' '
     out.append(namstr+'\n')
     out.append(valstr+'\n')
     out.append(sigstr+'\n')
-    return out
-
-def WriteSpnRBModel(RBModel,sigDict={},irb=None):
-    '''Write description of a spinning rigid body. Code here to make usable from G2export_CIF
-    '''
-    out = ' Radius: %10.3f refine? %s\n'%(RBModel['radius'][0],RBModel['radius'][1])    
     return out
 

trunk/GSASIIstrMain.py ¶

@@ -423 +423 @@
                 G2stIO.SetISOmodes(parmDict,sigDict,Phases,printFile)
                 G2stIO.SetHistogramPhaseData(parmDict,sigDict,Phases,Histograms,calcControls,
-                                             pFile=printFile,covMatrix=covMatrix,varyList=varyList)
+                    pFile=printFile,covMatrix=covMatrix,varyList=varyList)
                 G2stIO.SetHistogramData(parmDict,sigDict,Histograms,calcControls,pFile=printFile)
                 # check for variables outside their allowed range, reset and freeze them

trunk/GSASIIstrMath.py ¶

@@ -76 +76 @@
                     name = '::RBV;'+str(j)+':'+str(i)
                     VRBData[rbId]['VectMag'][j] = parmDict[name]
-    if RBIds['Spin']:
-        SRBData = RBData['Spin']
-        for i,rbId in enumerate(SRBIds):
-            if SRBData[rbId]['useCount']:
-                name = '::RBS;0:'+str(i)
-                SRBData[rbId]['radius'][0] = parmDict[name]
         
     for phase in Phases:
@@ -171 +165 @@
             jrb = SRBIds.index(RBObj['RBId'][0])
             name = pfx+'RBSOa:%d:%d'%(iAt,jrb)
-            RBObj['Orient'][0][0] = parmDict[name]
+            for i,po in enumerate(['RBSOa:','RBSOi:','RBSOj:','RBSOk:']):
+                name = pfx+'%s%d:%d'%(po,iAt,jrb)
+                RBObj['Orient'][0][i] = parmDict[name]                
             for ish in range(len(RBObj['RBId'])):
                 jrb = SRBIds.index(RBObj['RBId'][ish])
@@ -205 +201 @@
     RBModels =  Phase['RBModels']
     
-    for irb,RBObj in enumerate(RBModels.get('Spin',[])):
-        name = '::RBS;0:%d'%irb
-        if name in dFdvDict:
-            dFdvDict[name] += dFdvDict[pfx+'RBS;0:%d'%irb]
-        else:            
-            dFdvDict[name] = dFdvDict[pfx+'RBS;0:%d'%irb]
+    # for irb,RBObj in enumerate(RBModels.get('Spin',[])):
+    #     name = '::RBS;0:%d'%irb
+    #     if name in dFdvDict:
+    #         dFdvDict[name] += dFdvDict[pfx+'RBS;0:%d'%irb]
+    #     else:            
+    #         dFdvDict[name] = dFdvDict[pfx+'RBS;0:%d'%irb]
                    
     for irb,RBObj in enumerate(RBModels.get('Vector',[])):
@@ -399 +395 @@
                     break
                 Shell = SHdat[shl]
-                R = Shell['R']
+                R = Shell['Radius']
                 Atm = Shell['AtType']
                 Nat = Shell['Natoms']
@@ -409 +405 @@
                     dat = G2el.getBLvalues(BLtables)
                     SFF = dat[Atm]
-                Rname = ';0:%s'%Irb
+                Rname = 'Sh;%s;Radius:%d:%s'%(shl,iAt,Irb)
                 R0 = sp.spherical_jn(0,QR*R)/(4.*np.pi)
                 R0P = sp.spherical_jn(0,QR*(R+0.01))/(4.*np.pi)
@@ -440 +436 @@
                         name = 'Sh;%s;%s:%d:%s'%(shl,item,iAt,Irb)
                         dFFdS[name] = Nat*SFF*BS*SH
-                        #fill derivatives here wrt iAt,ishl,item or l,m
+                dFFdS[Rname] = dBSdR
                 ishl += 1
-                dFFdS[Rname] = dBSdR
             dFFdS[Oname] = dSHdO
             dFFdS[Oiname] = dSHdOi
@@ -468 +463 @@
             shno = bits[1]
             if shno not in SHCdict[atid]:
-                SHCdict[atid][shno] = {'R':parmDict['::RBS;0:%s'%shno]}
+                SHCdict[atid][shno] = {}
             if 'AtType' in bits[0] or 'Natoms' in bits[0] or 'ShR' in bits[0]:
                 SHCdict[atid][shno][bits[0]] = parmDict[parm]
@@ -1067 +1062 @@
     Tdata,Mdata,Fdata,Xdata,dXdata,IAdata,Uisodata,Uijdata,Gdata = \
         GetAtomFXU(pfx,calcControls,parmDict)
-    atFlg = np.zeros(len(Tdata))
+    atFlg = np.zeros(len(Tdata)) #non zero for Q type atoms - see below
     if not Xdata.size:          #no atoms in phase!
         return {}