Changeset 4425


Ignore:
Timestamp:
May 17, 2020 4:31:21 PM (17 months ago)
Author:
vondreele
Message:

fix generation of fullrmc pdb file for structures with partial substitutions & vacancies.

Location:
trunk
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIImath.py

    r4421 r4425  
    732732                IndB = ma.nonzero(ma.masked_greater(dist-radiusFactor*sumR[:,iA],0.))
    733733                for iU in IndB[0]:
    734                     if AtNames[iA] != notName:
     734                    if AtNames[iA%len(AtNames)] != notName:
    735735                        unit = Units[iU]
    736736                        if np.any(unit):
     
    739739                            Topstr = ' +(%4d)'%(Top)
    740740                        if Short:
    741                             Neigh.append([AtNames[iA],dist[iU],True])
     741                            Neigh.append([AtNames[iA%len(AtNames)],dist[iU],True])
    742742                        else:
    743743                            Neigh.append([AtNames[iA]+Topstr,atTypes[iA],dist[iU],dx[iU]])
  • trunk/GSASIIphsGUI.py

    r4424 r4425  
    55515551                                    NamesT = []
    55525552                                    for item in rdfDict:
     5553                                        if 'va' in item:
     5554                                            continue
    55535555                                        if 'rdf' not in item and 'g(r)' in title:
    55545556                                            Ylab = r'$\mathsf{G(r),\AA^{-1}}$'
  • trunk/GSASIIpwd.py

    r4420 r4425  
    22072207    return fname
    22082208
    2209 def MakeRMC6f(PWDdata,Name,Phase,RMCPdict):
    2210    
    2211     def findDup(Atoms):
    2212         Dup = []
    2213         Fracs = []
    2214         for iat1,at1 in enumerate(Atoms):
    2215             if any([at1[0] in dup for dup in Dup]):
    2216                 continue
    2217             else:
    2218                 Dup.append([at1[0],])
    2219                 Fracs.append([at1[6],])
    2220             for iat2,at2 in enumerate(Atoms[(iat1+1):]):
    2221                 if np.sum((np.array(at1[3:6])-np.array(at2[3:6]))**2) < 0.00001:
    2222                     Dup[-1] += [at2[0],]
    2223                     Fracs[-1]+= [at2[6],]
    2224         return Dup,Fracs
     2209def findDup(Atoms):
     2210    Dup = []
     2211    Fracs = []
     2212    for iat1,at1 in enumerate(Atoms):
     2213        if any([at1[0] in dup for dup in Dup]):
     2214            continue
     2215        else:
     2216            Dup.append([at1[0],])
     2217            Fracs.append([at1[6],])
     2218        for iat2,at2 in enumerate(Atoms[(iat1+1):]):
     2219            if np.sum((np.array(at1[3:6])-np.array(at2[3:6]))**2) < 0.00001:
     2220                Dup[-1] += [at2[0],]
     2221                Fracs[-1]+= [at2[6],]
     2222    return Dup,Fracs
     2223
     2224def MakeRMC6f(PWDdata,Name,Phase,RMCPdict):   
    22252225   
    22262226    Meta = RMCPdict['metadata']
     
    27532753    generalData = Phase['General']
    27542754    Atseq = RMCPdict['atSeq']
     2755    Dups,Fracs = findDup(Phase['Atoms'])
     2756    Sfracs = [np.cumsum(fracs) for fracs in Fracs]
    27552757    Supercell = RMCPdict['SuperCell']
    27562758    Cell = generalData['Cell'][1:7]
     
    27612763    newPhase,Atcodes = G2lat.TransformPhase(Phase,newPhase,Trans,np.zeros(3),np.zeros(3),ifMag=False,Force=False)
    27622764    Atoms = newPhase['Atoms']
    2763     XYZ = np.array([atom[3:6] for atom in Atoms]).T
     2765
     2766    Natm = np.core.defchararray.count(np.array(Atcodes),'+')    #no. atoms in original unit cell
     2767    Natm = np.count_nonzero(Natm-1)
     2768    Atoms = newPhase['Atoms']
     2769    Satoms = G2mth.sortArray(G2mth.sortArray(G2mth.sortArray(Atoms,5),4),3)
     2770    Datoms = [[atom for atom in Satoms if atom[0] in dup] for dup in Dups]
     2771    Natoms = []
     2772    for idup,dup in enumerate(Dups):
     2773        ldup = len(dup)
     2774        datoms = Datoms[idup]
     2775        natm = len(datoms)
     2776        i = 0
     2777        while i < natm:
     2778            atoms = datoms[i:i+ldup]
     2779            try:
     2780                atom = atoms[np.searchsorted(Sfracs[idup],rand.random())]
     2781                Natoms.append(atom)
     2782            except IndexError:      #what about vacancies?
     2783                if 'Va' not in Atseq:
     2784                    Atseq.append('Va')
     2785                    RMCPdict['aTypes']['Va'] = 0.0
     2786                atom = atoms[0]
     2787                atom[1] = 'Va'
     2788                Natoms.append(atom)
     2789            i += ldup
     2790
     2791
     2792
     2793
     2794
     2795    XYZ = np.array([atom[3:6] for atom in Natoms]).T
    27642796    XYZptp = np.array([ma.ptp(XYZ[0]),ma.ptp(XYZ[1]),ma.ptp(XYZ[2])])/2.
    27652797    Cell = newPhase['General']['Cell'][1:7]
     
    27802812    if RMCPdict['byMolec']:
    27812813        NPM = RMCPdict['Natoms']
    2782         for iat,atom in enumerate(Atoms):
     2814        for iat,atom in enumerate(Natoms):
    27832815            XYZ = np.inner(A,np.array(atom[3:6])-XYZptp)    #shift origin to middle & make Cartesian;residue = 'RMC'
    27842816            fl.write('ATOM  %5d %-4s RMC%6d%12.3f%8.3f%8.3f  1.00  0.00          %2s\n'%(       
     
    27872819    else:
    27882820        for atm in Atseq:
    2789             for iat,atom in enumerate(Atoms):
     2821            for iat,atom in enumerate(Natoms):
    27902822                if atom[1] == atm:
    27912823                    XYZ = np.inner(A,np.array(atom[3:6])-XYZptp)    #shift origin to middle & make Cartesian
  • trunk/imports/G2pdf_gr.py

    r4339 r4425  
    2323GSASIIpath.SetVersionNumber("$Revision$")
    2424
     25class txt_FSQReaderClass(G2obj.ImportPDFData):
     26    'Routines to import S(Q) data from a .fq file'
     27    def __init__(self):
     28        super(self.__class__,self).__init__( # fancy way to self-reference
     29            extensionlist=('.fq','.sq'),
     30            strictExtension=False,
     31            formatName = 'q (A-1) step S(Q) data',
     32            longFormatName = 'q (A-1) stepped S(Q) PDF data from pdfGet or GSAS-II'
     33            )
     34
     35    # Validate the contents -- make sure we only have valid lines
     36    def ContentsValidator(self, filename):
     37        'Look through the file for expected types of lines in a valid r-step file'
     38        filepointer = open(filename,'r')
     39        Ndata = 0
     40        for i,S in enumerate(filepointer):
     41            if '#' in S[:2]:
     42                break
     43            if len(S.split()) != 2:
     44                break
     45        for i,S in enumerate(filepointer):           
     46            vals = S.split()
     47            if len(vals) >= 2:
     48                try:
     49                    data = [float(val) for val in vals]
     50                    Ndata += 1
     51                except ValueError:
     52                    pass
     53        if not Ndata:     
     54            self.errors = 'No 2 or more column numeric data found'
     55            filepointer.close()
     56            return False
     57        filepointer.close()
     58        return True # no errors encountered
     59
     60    def Reader(self,filename,ParentFrame=None, **unused):
     61        print ('Read a q-step text file')
     62        x = []
     63        y = []
     64        ifData = False
     65        filepointer = open(filename,'r')
     66        for i,S in enumerate(filepointer):
     67            if not ifData:
     68                if len(S) == 1:     #skip blank line
     69                    continue
     70                if len(S.split()) != 2:
     71                    continue
     72                self.comments.append(S[:-1])
     73                if '#' in S[:2]:
     74                    ifData = True
     75            else:
     76                vals = S.split()
     77                if len(vals) >= 2:
     78                    try:
     79                        data = [float(val) for val in vals]
     80                        x.append(float(data[0]))
     81                        y.append(float(data[1]))
     82                    except ValueError:
     83                        msg = 'Error in line '+str(i+1)
     84                        print (msg)
     85                        continue
     86        self.pdfdata = np.array([
     87            np.array(x), # x-axis values q
     88            np.array(y), # pdf f(q))
     89            ])
     90        self.pdfentry[0] = filename
     91        self.pdfentry[2] = 1 # xy file only has one bank
     92        self.idstring = ospath.basename(filename)
     93
     94        return True
     95
    2596class txt_PDFReaderClass(G2obj.ImportPDFData):
    2697    'Routines to import PDF G(R) data from a .gr file'
     
    39110        Ndata = 0
    40111        for i,S in enumerate(filepointer):
    41             if '#L' in S[:2]:
     112            if '#' in S[:2]:
    42113                break
    43114        for i,S in enumerate(filepointer):           
     
    67138                    continue
    68139                self.comments.append(S[:-1])
    69                 if '#L' in S[:2]:
     140                if '#' in S[:2]:
    70141                    ifData = True
    71142            else:
Note: See TracChangeset for help on using the changeset viewer.