Changeset 943

Timestamp:

Jun 5, 2013 2:38:00 PM (10 years ago)

Author:

vondreele

Message:

more work on MC/SA

Location:

trunk

Files:

• GSASIImath.py (modified) (8 diffs)
• ReadMarCCDFrame.py (modified) (1 diff)

trunk/GSASIImath.py

@@ -1928 +1928 @@
         parmDict[pfx+'Amul'] = len(G2spc.GenAtom(XYZ,SGData))
             
-    def getRBparms(item,mfx,aTypes,RBdata,parmDict,varyList):
+    def getRBparms(item,mfx,aTypes,RBdata,atNo,parmDict,varyList):
         parmDict[mfx+'MolCent'] = item['MolCent']
         parmDict[mfx+'RBId'] = item['RBId']
@@ -1963 +1963 @@
                     lower.append(limits[0])
                     upper.append(limits[1])
-        aTypes |= set(RBdata[item['Type']][item['RBId']]['rbTypes'])
-
-    def GetAtomFX(pfx,Natoms,parmDict):
+        atypes = RBdata[item['Type']][item['RBId']]['rbTypes']
+        aTypes |= set(atypes)
+        atNo += len(atypes)
+        return atNo
+
+    def GetAtomTMX(pfx,RBdata,parmDict):
         'Needs a doc string'
-        Tdata = Natoms*[' ',]
-        Mdata = np.zeros(Natoms)
-        Fdata = np.zeros(Natoms)
-        Xdata = np.zeros((3,Natoms))
+        atNo = parmDIct['atNo']
+        nfixAt = parmDict['nfixAt']
+        Tdata = atNo*[' ',]
+        Mdata = np.zeros(atNo)
+        Fdata = np.zeros(atNo)
+        Xdata = np.zeros((3,atNo))
         keys = {'Atype:':Tdata,'Amul:':Mdata,
             'Ax:':Xdata[0],'Ay:':Xdata[1],'Az:':Xdata[2]}
-        for iatm in range(Natoms):
+        nObjs = parmDict['nObjs']
+        for iatm in range(nfixAt):
             for key in keys:
                 parm = pfx+key+str(iatm)
@@ -1980 +1986 @@
         return Tdata,Mdata,Xdata
     
-#    def mcsaSfCalc(refList,G,SGData,parmDict):
+#    def mcsaSfCalc(refList,RBdata,G,SGData,parmDict):
 #        ''' Compute structure factors for all h,k,l for phase
 #        input:
@@ -2036 +2042 @@
     generalData = data['General']
     SGData = generalData['SGData']
-    pId = data['pId']
     fixAtoms = data['Atoms']                       #if any
     cx,ct,cs = generalData['AtomPtrs'][:3]
@@ -2052 +2057 @@
             name = pfx+pstr[i]
             parmDict[name] = atm[cx+i]
-        atNo += 1        
+        atNo += 1
+    parmDict['nfixAt'] = len(fixAtoms)        
     mcsaControls = generalData['MCSA controls']
     reflName = mcsaControls['Data source']
@@ -2061 +2067 @@
     for i,item in enumerate(MCSAObjs):
         mfx = str(i)+':'
+        parmDict[mfx+'Type'] = item['Type']
         if item['Type'] == 'MD':
             getMDparms(item,mfx,parmDict,varyList)
@@ -2069 +2076 @@
         elif item['Type'] in ['Residue','Vector']:
             pfx = mfx+':'
-            getRBparms(item,pfx,aTypes,RBdata,parmDict,varyList)
+            atNo = getRBparms(item,pfx,aTypes,RBdata,atNo,parmDict,varyList)
+    parmDict['atNo'] = atNo                 #total no. of atoms
+    parmDict['nObj'] = len(MCSAObjs)
     FFtables = G2el.GetFFtable(aTypes)
     refs = []
     if 'PWDR' in reflName:
         for ref in reflData:
-            h,k,l,m,d,pos,sig,gam,Fsq = ref[:9]
+            h,k,l,m,d,pos,sig,gam,f = ref[:9]
             if d >= mcsaControls['dmin']:
                 sig = gamFW(sig,gam)/sq8ln2        #--> sig from FWHM
@@ -2105 +2114 @@
                 Uniq,phi = G2spc.GenHKLf([h,k,l],SGData)[2:]
                 FFs = G2el.getFFvalues(FFtables,SQ)
-                refs.append([h,k,l,m,f*m,0.,0.,FFs,Uniq,phi])
+                refs.append([h,k,l,m,f*m,iref,0.,FFs,Uniq,phi])
         nRef = len(refs)
+        pfx = str(data['pId'])+'::PWLref:'
         rcov = np.zeros((nRef,nRef))        
-        Iref = 0
-        for iref,refI in enumerate(reflData):
-            if refI[4] >= mcsaControls['dmin'] and refI[5]:       #skip unrefined ones
-                nameI = pfx+'PWLref:'+str(iref)
-                if nameI in covData['varyList']:
-                    Iindx = vList.index(nameI)
-                    rcov[Iref][Iref] = covMatrix[Iindx][Iindx]
-                    Jref = 0
-                    for jref,refJ in enumerate(reflData[:iref]):
-                        if refJ[5]:
-                            nameJ = pfx+'PWLref:'+str(jref)
-                            try:
-                                rcov[Iref][Jref] = covMatrix[vList.index(nameI)][vList.index(nameJ)]
-                            except ValueError:
-                                rcov[Iref][Jref] = rcov[Iref][Jref-1]
-                            Jref += 1
-                else:
-                    rcov[Iref] = rcov[Iref-1]
-                    rcov[Iref][Iref] = rcov[Iref-1][Iref-1]
-                Iref += 1
+        for iref,refI in enumerate(refs):
+            I = refI[5]
+            nameI = pfx+str(I)
+            if nameI in vList:
+                Iindx = vList.index(nameI)
+                rcov[iref][iref] = covMatrix[Iindx][Iindx]
+                for jref,refJ in enumerate(refs[:iref]):
+                    J = refJ[5]
+                    nameJ = pfx+str(J)
+                    try:
+                        rcov[iref][jref] = covMatrix[vList.index(nameI)][vList.index(nameJ)]
+                    except ValueError:
+                        rcov[iref][jref] = rcov[iref][jref-1]
+            else:
+                rcov[iref] = rcov[iref-1]
+                rcov[iref][iref] = rcov[iref-1][iref-1]
         rcov += (rcov.T-np.diagflat(np.diagonal(rcov)))
         Rdiag = np.sqrt(np.diag(rcov))

trunk/ReadMarCCDFrame.py

@@ -446 +446 @@
         File.seek(4096)
         self.image = ar.array(self.arrayTypeCode,File.read(self.byteDepth*self.TIFFsizeX*self.TIFFsizeY))
-        # reverse the array so if can have the same view as is read in marccd
-        # also switch the view direction 
-        self.image.reverse()
+    # reverse the array so if can have the same view as is read in marccd
+    # also switch the view direction 
+        self.image.reverse()
         self.viewDirection = abs(self.viewDirection - 1)