Changeset 3990 for trunk/GSASIIpwd.py
- Timestamp:
- May 22, 2019 2:04:11 PM (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/GSASIIpwd.py
r3959 r3990 2810 2810 Layers['Sadp']['Img'] = Sapd 2811 2811 print (' GETSAD time = %.2fs'%(time.time()-time0)) 2812 # GSASIIpath.IPyBreak() 2813 2812 2813 ############################################################################### 2814 #### Maximum Entropy Method - Dysnomia 2815 ############################################################################### 2816 2817 def makePRFfile(data,MEMtype): 2818 ''' makes Dysnomia .prf control file from Dysnomia GUI controls 2819 2820 ;param dict data: GSAS-II phase data 2821 :param int MEMtype: 1 for neutron data with negative scattering lengths 2822 0 otherwise 2823 :returns str: name of Dysnomia control file 2824 ''' 2825 2826 generalData = data['General'] 2827 pName = generalData['Name'].replace(' ','_') 2828 DysData = data['Dysnomia'] 2829 prfName = pName+'.prf' 2830 prf = open(prfName,'w') 2831 prf.write('$PREFERENCES\n') 2832 prf.write(pName+'.mem\n') #or .fos? 2833 prf.write(pName+'.out\n') 2834 prf.write(pName+'.pgrid\n') 2835 prf.write(pName+'.fba\n') 2836 prf.write(pName+'_eps.raw\n') 2837 prf.write('%d\n'%MEMtype) 2838 if DysData['DenStart'] == 'uniform': 2839 prf.write('0\n') 2840 else: 2841 prf.write('1\n') 2842 if DysData['Optimize'] == 'ZSPA': 2843 prf.write('0\n') 2844 else: 2845 prf.write('1\n') 2846 prf.write('1\n') 2847 if DysData['Lagrange'][0] == 'user': 2848 prf.write('0\n') 2849 else: 2850 prf.write('1\n') 2851 prf.write('%.4f %d\n'%(DysData['Lagrange'][1],DysData['wt pwr'])) 2852 prf.write('%.3f\n'%DysData['Lagrange'][2]) 2853 prf.write('%.2f\n'%DysData['E_factor']) 2854 prf.write('1\n') 2855 prf.write('0\n') 2856 prf.write('%d\n'%DysData['Ncyc']) 2857 prf.write('1\n') 2858 prf.write('1 0 0 0 0 0 0 0\n') 2859 if DysData['prior'] == 'uniform': 2860 prf.write('0\n') 2861 else: 2862 prf.write('1\n') 2863 prf.close() 2864 return prfName 2865 2866 def makeMEMfile(data,reflData,MEMtype): 2867 ''' make Dysnomia .mem file of reflection data, etc. 2868 ;param dict data: GSAS-II phase data 2869 :param list reflData: GSAS-II reflection data 2870 :param int MEMtype: 1 for neutron data with negative scattering lengths 2871 0 otherwise 2872 ''' 2873 2874 DysData = data['Dysnomia'] 2875 generalData = data['General'] 2876 pName = generalData['Name'].replace(' ','_') 2877 memName = pName+'.mem' 2878 Map = generalData['Map'] 2879 Type = Map['Type'] 2880 UseList = Map['RefList'] 2881 mem = open(memName,'w') 2882 mem.write('%s\n'%(generalData['Name']+' from '+UseList[0])) 2883 a,b,c,alp,bet,gam = generalData['Cell'][1:7] 2884 mem.write('%10.5f%10.5f%10.5f%10.5f,%10.5f%10.5f\n'%(a,b,c,alp,bet,gam)) 2885 mem.write(' 0.0000000 0.0000000 -1 0 0 0 P\n') #dummy PO stuff 2886 SGSym = generalData['SGData']['SpGrp'] 2887 try: 2888 SGId = G2spc.spgbyNum.index(SGSym) 2889 except IndexError: 2890 return False 2891 org = 1 2892 if SGSym in G2spc.spg2origins: 2893 org = 2 2894 mapsize = Map['rho'].shape 2895 sumZ = 0. 2896 sumpos = 0. 2897 sumneg = 0. 2898 mem.write('%5d%5d%5d%5d%5d\n'%(SGId,org,mapsize[0],mapsize[1],mapsize[2])) 2899 for atm in generalData['NoAtoms']: 2900 Nat = generalData['NoAtoms'][atm] 2901 AtInfo = G2elem.GetAtomInfo(atm) 2902 sumZ += Nat*AtInfo['Z'] 2903 isotope = generalData['Isotope'][atm] 2904 blen = generalData['Isotopes'][atm][isotope]['SL'][0] 2905 if blen < 0.: 2906 sumneg += blen*Nat 2907 else: 2908 sumpos += blen*Nat 2909 if 'X' in Type: 2910 mem.write('%10.2f 0.001\n'%sumZ) 2911 elif 'N' in Type and MEMtype: 2912 mem.write('%10.3f%10.3f 0.001\n'%(sumpos,sumneg)) 2913 else: 2914 mem.write('%10.3 0.001\n'%sumpos) 2915 2916 refs = [] 2917 prevpos = 0. 2918 for ref in reflData: 2919 if 'T' in Type: 2920 h,k,l,mult,dsp,pos,sig,gam,Fobs,Fcalc,phase,x,x,x,x,prfo = ref[:16] 2921 FWHM = getgamFW(gam,sig) 2922 else: 2923 h,k,l,mult,dsp,pos,sig,gam,Fobs,Fcalc,phase,x,prfo = ref[:13] 2924 g = gam/100. #centideg -> deg 2925 s = np.sqrt(max(sig,0.0001))/100. #var -> sig in deg 2926 FWHM = getgamFW(g,s) 2927 delt = pos-prevpos 2928 refs.append([h,k,l,mult,pos,FWHM,Fobs,phase,delt]) 2929 prevpos = pos 2930 2931 ovlp = DysData['overlap'] 2932 refs1 = [] 2933 refs2 = [] 2934 nref2 = 0 2935 iref = 0 2936 Nref = len(refs) 2937 start = False 2938 while iref < Nref-1: 2939 if refs[iref+1][-1] < ovlp*refs[iref][5]: 2940 if refs[iref][-1] > ovlp*refs[iref][5]: 2941 refs2.append([]) 2942 start = True 2943 if nref2 == len(refs2): 2944 refs2.append([]) 2945 refs2[nref2].append(refs[iref]) 2946 else: 2947 if start: 2948 refs2[nref2].append(refs[iref]) 2949 start = False 2950 nref2 += 1 2951 else: 2952 refs1.append(refs[iref]) 2953 iref += 1 2954 if start: 2955 refs2[nref2].append(refs[iref]) 2956 else: 2957 refs1.append(refs[iref]) 2958 2959 mem.write('%5d\n'%len(refs1)) 2960 for ref in refs1: 2961 h,k,l = ref[:3] 2962 Fobs = np.sqrt(ref[6]) 2963 mem.write('%5d%5d%5d%10.3f%10.3f%10.3f\n'%(h,k,l,Fobs*npcosd(ref[7]),Fobs*npsind(ref[7]),max(0.01*Fobs,0.1))) 2964 mem.write('%5d\n'%len(refs2)) 2965 for ref2 in refs2: 2966 if not len(ref2): 2967 break 2968 mem.write('%5d\n'%len(ref2)) 2969 Gsum = 0. 2970 Msum = 0 2971 for ref in ref2: 2972 Gsum += ref[6]*ref[3] 2973 Msum += ref[3] 2974 G = np.sqrt(Gsum/Msum) 2975 h,k,l = ref2[0][:3] 2976 mem.write('%5d%5d%5d%10.3f%10.3f%5d\n'%(h,k,l,G,max(0.01*G,0.1),ref2[0][3])) 2977 for ref in ref2[1:]: 2978 h,k,l,m = ref[:4] 2979 mem.write('%5d%5d%5d%5d\n'%(h,k,l,m)) 2980 mem.write('0\n') 2981 mem.close() 2982 return True 2814 2983 #testing data 2815 2984 NeedTestData = True
Note: See TracChangeset
for help on using the changeset viewer.