Context Navigation

source: trunk/GSASIIstrIO.py @ 1570

Last change on this file since 1570 was 1570, checked in by vondreele, 8 years ago
fix to G2strIO for space group output. finish interpretation of super symmetry symbols
Property svn:eol-style set to `native` Property svn:keywords set to `Date Author Revision URL Id`
File size: 114.1 KB

Line
1	# -- coding: utf-8 --
2	########### SVN repository information ###################
3	# $Date: 2014-11-13 21:36:54 +0000 (Thu, 13 Nov 2014) $
4	# $Author: vondreele $
5	# $Revision: 1570 $
6	# $URL: trunk/GSASIIstrIO.py $
7	# $Id: GSASIIstrIO.py 1570 2014-11-13 21:36:54Z vondreele $
8	########### SVN repository information ###################
9	'''
10	GSASIIstrIO: structure I/O routines
11	-------------------------------------
12
13	'''
14	import sys
15	import os
16	import os.path as ospath
17	import time
18	import math
19	import copy
20	import random
21	import cPickle
22	import numpy as np
23	import numpy.ma as ma
24	import numpy.linalg as nl
25	import scipy.optimize as so
26	import GSASIIpath
27	GSASIIpath.SetVersionNumber("$Revision: 1570 $")
28	import GSASIIElem as G2el
29	import GSASIIgrid as G2gd
30	import GSASIIlattice as G2lat
31	import GSASIIspc as G2spc
32	import GSASIIobj as G2obj
33	import GSASIImapvars as G2mv
34	import GSASIImath as G2mth
35
36	sind = lambda x: np.sin(x*np.pi/180.)
37	cosd = lambda x: np.cos(x*np.pi/180.)
38	tand = lambda x: np.tan(x*np.pi/180.)
39	asind = lambda x: 180.*np.arcsin(x)/np.pi
40	acosd = lambda x: 180.*np.arccos(x)/np.pi
41	atan2d = lambda y,x: 180.*np.arctan2(y,x)/np.pi
42
43	ateln2 = 8.0*math.log(2.0)
44
45	def GetControls(GPXfile):
46	''' Returns dictionary of control items found in GSASII gpx file
47
48	:param str GPXfile: full .gpx file name
49	:return: dictionary of control items
50	'''
51	Controls = copy.copy(G2obj.DefaultControls)
52	fl = open(GPXfile,'rb')
53	while True:
54	try:
55	data = cPickle.load(fl)
56	except EOFError:
57	break
58	datum = data[0]
59	if datum[0] == 'Controls':
60	Controls.update(datum[1])
61	fl.close()
62	return Controls
63
64	def GetConstraints(GPXfile):
65	'''Read the constraints from the GPX file and interpret them
66
67	called in :func:`CheckConstraints`, :func:`GSASIIstrMain.Refine`
68	and :func:`GSASIIstrMain.SeqRefine`.
69	'''
70	fl = open(GPXfile,'rb')
71	while True:
72	try:
73	data = cPickle.load(fl)
74	except EOFError:
75	break
76	datum = data[0]
77	if datum[0] == 'Constraints':
78	constList = []
79	for item in datum[1]:
80	if item.startswith('_'): continue
81	constList += datum[1][item]
82	fl.close()
83	constDict,fixedList,ignored = ProcessConstraints(constList)
84	if ignored:
85	print ignored,'Constraints were rejected. Was a constrained phase, histogram or atom deleted?'
86	return constDict,fixedList
87	fl.close()
88	raise Exception,"No constraints in GPX file"
89
90	def ProcessConstraints(constList):
91	"""Interpret the constraints in the constList input into a dictionary, etc.
92	All :class:`GSASIIobj.G2VarObj` objects are mapped to the appropriate
93	phase/hist/atoms based on the object internals (random Ids). If this can't be
94	done (if a phase has been deleted, etc.), the variable is ignored.
95	If the constraint cannot be used due to too many dropped variables,
96	it is counted as ignored.
97
98	:param list constList: a list of lists where each item in the outer list
99	specifies a constraint of some form, as described in the :mod:`GSASIIobj`
100	:ref:`Constraint definition <Constraint_definitions_table>`.
101
102	:returns: a tuple of (constDict,fixedList,ignored) where:
103
104	* constDict (list of dicts) contains the constraint relationships
105	* fixedList (list) contains the fixed values for each type
106	of constraint.
107	* ignored (int) counts the number of invalid constraint items
108	(should always be zero!)
109	"""
110	constDict = []
111	fixedList = []
112	ignored = 0
113	for item in constList:
114	if item[-1] == 'h':
115	# process a hold
116	fixedList.append('0')
117	var = str(item[0][1])
118	if '?' not in var:
119	constDict.append({var:0.0})
120	else:
121	ignored += 1
122	elif item[-1] == 'f':
123	# process a new variable
124	fixedList.append(None)
125	D = {}
126	varyFlag = item[-2]
127	varname = item[-3]
128	for term in item[:-3]:
129	var = str(term[1])
130	if '?' not in var:
131	D[var] = term[0]
132	if len(D) > 1:
133	# add extra dict terms for input variable name and vary flag
134	if varname is not None:
135	if varname.startswith('::'):
136	varname = varname[2:].replace(':',';')
137	else:
138	varname = varname.replace(':',';')
139	D['_name'] = '::' + varname
140	D['_vary'] = varyFlag == True # force to bool
141	constDict.append(D)
142	else:
143	ignored += 1
144	#constFlag[-1] = ['Vary']
145	elif item[-1] == 'c':
146	# process a contraint relationship
147	D = {}
148	for term in item[:-3]:
149	var = str(term[1])
150	if '?' not in var:
151	D[var] = term[0]
152	if len(D) >= 1:
153	fixedList.append(str(item[-3]))
154	constDict.append(D)
155	else:
156	ignored += 1
157	elif item[-1] == 'e':
158	# process an equivalence
159	firstmult = None
160	eqlist = []
161	for term in item[:-3]:
162	if term[0] == 0: term[0] = 1.0
163	var = str(term[1])
164	if '?' in var: continue
165	if firstmult is None:
166	firstmult = term[0]
167	firstvar = var
168	else:
169	eqlist.append([var,firstmult/term[0]])
170	if len(eqlist) > 0:
171	G2mv.StoreEquivalence(firstvar,eqlist)
172	else:
173	ignored += 1
174	else:
175	ignored += 1
176	return constDict,fixedList,ignored
177
178	def CheckConstraints(GPXfile):
179	'''Load constraints and related info and return any error or warning messages'''
180	# init constraints
181	G2mv.InitVars()
182	# get variables
183	Histograms,Phases = GetUsedHistogramsAndPhases(GPXfile)
184	if not Phases:
185	return 'Error: No Phases!',''
186	if not Histograms:
187	return 'Error: no diffraction data',''
188	rigidbodyDict = GetRigidBodies(GPXfile)
189	rbIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[]})
190	rbVary,rbDict = GetRigidBodyModels(rigidbodyDict,Print=False)
191	Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables = GetPhaseData(Phases,RestraintDict=None,rbIds=rbIds,Print=False)
192	hapVary,hapDict,controlDict = GetHistogramPhaseData(Phases,Histograms,Print=False)
193	histVary,histDict,controlDict = GetHistogramData(Histograms,Print=False)
194	varyList = rbVary+phaseVary+hapVary+histVary
195	constrDict,fixedList = GetConstraints(GPXfile)
196	return G2mv.CheckConstraints(varyList,constrDict,fixedList)
197
198	def GetRestraints(GPXfile):
199	'''Read the restraints from the GPX file.
200	Throws an exception if not found in the .GPX file
201	'''
202	fl = open(GPXfile,'rb')
203	while True:
204	try:
205	data = cPickle.load(fl)
206	except EOFError:
207	break
208	datum = data[0]
209	if datum[0] == 'Restraints':
210	restraintDict = datum[1]
211	fl.close()
212	return restraintDict
213
214	def GetRigidBodies(GPXfile):
215	'''Read the rigid body models from the GPX file
216	'''
217	fl = open(GPXfile,'rb')
218	while True:
219	try:
220	data = cPickle.load(fl)
221	except EOFError:
222	break
223	datum = data[0]
224	if datum[0] == 'Rigid bodies':
225	rigidbodyDict = datum[1]
226	fl.close()
227	return rigidbodyDict
228
229	def GetFprime(controlDict,Histograms):
230	'Needs a doc string'
231	FFtables = controlDict['FFtables']
232	if not FFtables:
233	return
234	histoList = Histograms.keys()
235	histoList.sort()
236	for histogram in histoList:
237	if histogram[:4] in ['PWDR','HKLF']:
238	Histogram = Histograms[histogram]
239	hId = Histogram['hId']
240	hfx = ':%d:'%(hId)
241	if 'X' in controlDict[hfx+'histType']:
242	keV = controlDict[hfx+'keV']
243	for El in FFtables:
244	Orbs = G2el.GetXsectionCoeff(El.split('+')[0].split('-')[0])
245	FP,FPP,Mu = G2el.FPcalc(Orbs, keV)
246	FFtables[El][hfx+'FP'] = FP
247	FFtables[El][hfx+'FPP'] = FPP
248
249	def GetPhaseNames(GPXfile):
250	''' Returns a list of phase names found under 'Phases' in GSASII gpx file
251
252	:param str GPXfile: full .gpx file name
253	:return: list of phase names
254	'''
255	fl = open(GPXfile,'rb')
256	PhaseNames = []
257	while True:
258	try:
259	data = cPickle.load(fl)
260	except EOFError:
261	break
262	datum = data[0]
263	if 'Phases' == datum[0]:
264	for datus in data[1:]:
265	PhaseNames.append(datus[0])
266	fl.close()
267	return PhaseNames
268
269	def GetAllPhaseData(GPXfile,PhaseName):
270	''' Returns the entire dictionary for PhaseName from GSASII gpx file
271
272	:param str GPXfile: full .gpx file name
273	:param str PhaseName: phase name
274	:return: phase dictionary
275	'''
276	fl = open(GPXfile,'rb')
277	General = {}
278	Atoms = []
279	while True:
280	try:
281	data = cPickle.load(fl)
282	except EOFError:
283	break
284	datum = data[0]
285	if 'Phases' == datum[0]:
286	for datus in data[1:]:
287	if datus[0] == PhaseName:
288	break
289	fl.close()
290	return datus[1]
291
292	def GetHistograms(GPXfile,hNames):
293	""" Returns a dictionary of histograms found in GSASII gpx file
294
295	:param str GPXfile: full .gpx file name
296	:param str hNames: list of histogram names
297	:return: dictionary of histograms (types = PWDR & HKLF)
298
299	"""
300	fl = open(GPXfile,'rb')
301	Histograms = {}
302	while True:
303	try:
304	data = cPickle.load(fl)
305	except EOFError:
306	break
307	datum = data[0]
308	hist = datum[0]
309	if hist in hNames:
310	if 'PWDR' in hist[:4]:
311	PWDRdata = {}
312	PWDRdata.update(datum[1][0]) #weight factor
313	PWDRdata['Data'] = ma.array(ma.getdata(datum[1][1])) #masked powder data arrays/clear previous masks
314	PWDRdata[data[2][0]] = data[2][1] #Limits & excluded regions (if any)
315	PWDRdata[data[3][0]] = data[3][1] #Background
316	PWDRdata[data[4][0]] = data[4][1] #Instrument parameters
317	PWDRdata[data[5][0]] = data[5][1] #Sample parameters
318	try:
319	PWDRdata[data[9][0]] = data[9][1] #Reflection lists might be missing
320	except IndexError:
321	PWDRdata['Reflection Lists'] = {}
322	PWDRdata['Residuals'] = {}
323
324	Histograms[hist] = PWDRdata
325	elif 'HKLF' in hist[:4]:
326	HKLFdata = {}
327	HKLFdata.update(datum[1][0]) #weight factor
328	#patch
329	if 'list' in str(type(datum[1][1])):
330	#if isinstance(datum[1][1],list):
331	RefData = {'RefList':[],'FF':{}}
332	for ref in datum[1][1]:
333	RefData['RefList'].append(ref[:11]+[ref[13],])
334	RefData['RefList'] = np.array(RefData['RefList'])
335	datum[1][1] = RefData
336	#end patch
337	datum[1][1]['FF'] = {}
338	HKLFdata['Data'] = datum[1][1]
339	HKLFdata[data[1][0]] = data[1][1] #Instrument parameters
340	HKLFdata['Reflection Lists'] = None
341	HKLFdata['Residuals'] = {}
342	Histograms[hist] = HKLFdata
343	fl.close()
344	return Histograms
345
346	def GetHistogramNames(GPXfile,hType):
347	""" Returns a list of histogram names found in GSASII gpx file
348
349	:param str GPXfile: full .gpx file name
350	:param str hType: list of histogram types
351	:return: list of histogram names (types = PWDR & HKLF)
352
353	"""
354	fl = open(GPXfile,'rb')
355	HistogramNames = []
356	while True:
357	try:
358	data = cPickle.load(fl)
359	except EOFError:
360	break
361	datum = data[0]
362	if datum[0][:4] in hType:
363	HistogramNames.append(datum[0])
364	fl.close()
365	return HistogramNames
366
367	def GetUsedHistogramsAndPhases(GPXfile):
368	''' Returns all histograms that are found in any phase
369	and any phase that uses a histogram. This also
370	assigns numbers to used phases and histograms by the
371	order they appear in the file.
372
373	:param str GPXfile: full .gpx file name
374	:returns: (Histograms,Phases)
375
376	* Histograms = dictionary of histograms as {name:data,...}
377	* Phases = dictionary of phases that use histograms
378
379	'''
380	phaseNames = GetPhaseNames(GPXfile)
381	histoList = GetHistogramNames(GPXfile,['PWDR','HKLF'])
382	allHistograms = GetHistograms(GPXfile,histoList)
383	phaseData = {}
384	for name in phaseNames:
385	phaseData[name] = GetAllPhaseData(GPXfile,name)
386	Histograms = {}
387	Phases = {}
388	for phase in phaseData:
389	Phase = phaseData[phase]
390	if Phase['Histograms']:
391	if phase not in Phases:
392	pId = phaseNames.index(phase)
393	Phase['pId'] = pId
394	Phases[phase] = Phase
395	for hist in Phase['Histograms']:
396	if 'Use' not in Phase['Histograms'][hist]: #patch
397	Phase['Histograms'][hist]['Use'] = True
398	if hist not in Histograms and Phase['Histograms'][hist]['Use']:
399	try:
400	Histograms[hist] = allHistograms[hist]
401	hId = histoList.index(hist)
402	Histograms[hist]['hId'] = hId
403	except KeyError: # would happen if a referenced histogram were
404	# renamed or deleted
405	print('For phase "'+str(phase)+
406	'" unresolved reference to histogram "'+str(hist)+'"')
407	G2obj.IndexAllIds(Histograms=Histograms,Phases=Phases)
408	return Histograms,Phases
409
410	def getBackupName(GPXfile,makeBack):
411	'''
412	Get the name for the backup .gpx file name
413
414	:param str GPXfile: full .gpx file name
415	:param bool makeBack: if True the name of a new file is returned, if
416	False the name of the last file that exists is returned
417	:returns: the name of a backup file
418
419	'''
420	GPXpath,GPXname = ospath.split(GPXfile)
421	if GPXpath == '': GPXpath = '.'
422	Name = ospath.splitext(GPXname)[0]
423	files = os.listdir(GPXpath)
424	last = 0
425	for name in files:
426	name = name.split('.')
427	if len(name) == 3 and name[0] == Name and 'bak' in name[1]:
428	if makeBack:
429	last = max(last,int(name[1].strip('bak'))+1)
430	else:
431	last = max(last,int(name[1].strip('bak')))
432	GPXback = ospath.join(GPXpath,ospath.splitext(GPXname)[0]+'.bak'+str(last)+'.gpx')
433	return GPXback
434
435	def GPXBackup(GPXfile,makeBack=True):
436	'''
437	makes a backup of the current .gpx file (?)
438
439	:param str GPXfile: full .gpx file name
440	:param bool makeBack: if True (default), the backup is written to
441	a new file; if False, the last backup is overwritten
442	:returns: the name of the backup file that was written
443	'''
444	import distutils.file_util as dfu
445	GPXback = getBackupName(GPXfile,makeBack)
446	dfu.copy_file(GPXfile,GPXback)
447	return GPXback
448
449	def SetUsedHistogramsAndPhases(GPXfile,Histograms,Phases,RigidBodies,CovData,makeBack=True):
450	''' Updates gpxfile from all histograms that are found in any phase
451	and any phase that used a histogram. Also updates rigid body definitions.
452
453
454	:param str GPXfile: full .gpx file name
455	:param dict Histograms: dictionary of histograms as {name:data,...}
456	:param dict Phases: dictionary of phases that use histograms
457	:param dict RigidBodies: dictionary of rigid bodies
458	:param dict CovData: dictionary of refined variables, varyList, & covariance matrix
459	:param bool makeBack: True if new backup of .gpx file is to be made; else use the last one made
460
461	'''
462
463	GPXback = GPXBackup(GPXfile,makeBack)
464	print 'Read from file:',GPXback
465	print 'Save to file :',GPXfile
466	infile = open(GPXback,'rb')
467	outfile = open(GPXfile,'wb')
468	while True:
469	try:
470	data = cPickle.load(infile)
471	except EOFError:
472	break
473	datum = data[0]
474	# print 'read: ',datum[0]
475	if datum[0] == 'Phases':
476	for iphase in range(len(data)):
477	if data[iphase][0] in Phases:
478	phaseName = data[iphase][0]
479	data[iphase][1].update(Phases[phaseName])
480	elif datum[0] == 'Covariance':
481	data[0][1] = CovData
482	elif datum[0] == 'Rigid bodies':
483	data[0][1] = RigidBodies
484	try:
485	histogram = Histograms[datum[0]]
486	# print 'found ',datum[0]
487	data[0][1][1] = histogram['Data']
488	data[0][1][0].update(histogram['Residuals'])
489	for datus in data[1:]:
490	# print ' read: ',datus[0]
491	if datus[0] in ['Background','Instrument Parameters','Sample Parameters','Reflection Lists']:
492	datus[1] = histogram[datus[0]]
493	except KeyError:
494	pass
495
496	cPickle.dump(data,outfile,1)
497	infile.close()
498	outfile.close()
499	print 'GPX file save successful'
500
501	def SetSeqResult(GPXfile,Histograms,SeqResult):
502	'''
503	Needs doc string
504
505	:param str GPXfile: full .gpx file name
506	'''
507	GPXback = GPXBackup(GPXfile)
508	print 'Read from file:',GPXback
509	print 'Save to file :',GPXfile
510	infile = open(GPXback,'rb')
511	outfile = open(GPXfile,'wb')
512	while True:
513	try:
514	data = cPickle.load(infile)
515	except EOFError:
516	break
517	datum = data[0]
518	if datum[0] == 'Sequential results':
519	data[0][1] = SeqResult
520	# reset the Copy Next flag, since it should not be needed twice in a row
521	if datum[0] == 'Controls':
522	data[0][1]['Copy2Next'] = False
523	try:
524	histogram = Histograms[datum[0]]
525	data[0][1][1] = list(histogram['Data'])
526	for datus in data[1:]:
527	if datus[0] in ['Background','Instrument Parameters','Sample Parameters','Reflection Lists']:
528	datus[1] = histogram[datus[0]]
529	except KeyError:
530	pass
531
532	cPickle.dump(data,outfile,1)
533	infile.close()
534	outfile.close()
535	print 'GPX file save successful'
536
537	def ShowBanner(pFile=None):
538	'Print authorship, copyright and citation notice'
539	print >>pFile,80''
540	print >>pFile,' General Structure Analysis System-II Crystal Structure Refinement'
541	print >>pFile,' by Robert B. Von Dreele & Brian H. Toby'
542	print >>pFile,' Argonne National Laboratory(C), 2010'
543	print >>pFile,' This product includes software developed by the UChicago Argonne, LLC,'
544	print >>pFile,' as Operator of Argonne National Laboratory.'
545	print >>pFile,' Please cite:'
546	print >>pFile,' B.H. Toby & R.B. Von Dreele, J. Appl. Cryst. 46, 544-549 (2013)'
547
548	print >>pFile,80'','\n'
549
550	def ShowControls(Controls,pFile=None,SeqRef=False):
551	'Print controls information'
552	print >>pFile,' Least squares controls:'
553	print >>pFile,' Refinement type: ',Controls['deriv type']
554	if 'Hessian' in Controls['deriv type']:
555	print >>pFile,' Maximum number of cycles:',Controls['max cyc']
556	else:
557	print >>pFile,' Minimum delta-M/M for convergence: ','%.2g'%(Controls['min dM/M'])
558	print >>pFile,' Initial shift factor: ','%.3f'%(Controls['shift factor'])
559	if SeqRef:
560	print >>pFile,' Sequential refinement controls:'
561	print >>pFile,' Copy of histogram results to next: ',Controls['Copy2Next']
562	print >>pFile,' Process histograms in reverse order: ',Controls['Reverse Seq']
563
564	def GetPawleyConstr(SGLaue,PawleyRef,pawleyVary):
565	'needs a doc string'
566	# if SGLaue in ['-1','2/m','mmm']:
567	# return #no Pawley symmetry required constraints
568	eqvDict = {}
569	for i,varyI in enumerate(pawleyVary):
570	eqvDict[varyI] = []
571	refI = int(varyI.split(':')[-1])
572	ih,ik,il = PawleyRef[refI][:3]
573	dspI = PawleyRef[refI][4]
574	for varyJ in pawleyVary[i+1:]:
575	refJ = int(varyJ.split(':')[-1])
576	jh,jk,jl = PawleyRef[refJ][:3]
577	dspJ = PawleyRef[refJ][4]
578	if SGLaue in ['4/m','4/mmm']:
579	isum = ih2+ik2
580	jsum = jh2+jk2
581	if abs(il) == abs(jl) and isum == jsum:
582	eqvDict[varyI].append(varyJ)
583	elif SGLaue in ['3R','3mR']:
584	isum = ih2+ik2+il**2
585	jsum = jh2+jk2+jl**2
586	isum2 = ihik+ihil+ik*il
587	jsum2 = jhjk+jhjl+jk*jl
588	if isum == jsum and isum2 == jsum2:
589	eqvDict[varyI].append(varyJ)
590	elif SGLaue in ['3','3m1','31m','6/m','6/mmm']:
591	isum = ih2+ik2+ih*ik
592	jsum = jh2+jk2+jh*jk
593	if abs(il) == abs(jl) and isum == jsum:
594	eqvDict[varyI].append(varyJ)
595	elif SGLaue in ['m3','m3m']:
596	isum = ih2+ik2+il**2
597	jsum = jh2+jk2+jl**2
598	if isum == jsum:
599	eqvDict[varyI].append(varyJ)
600	elif abs(dspI-dspJ)/dspI < 1.e-4:
601	eqvDict[varyI].append(varyJ)
602	for item in pawleyVary:
603	if eqvDict[item]:
604	for item2 in pawleyVary:
605	if item2 in eqvDict[item]:
606	eqvDict[item2] = []
607	G2mv.StoreEquivalence(item,eqvDict[item])
608
609	def cellVary(pfx,SGData):
610	'needs a doc string'
611	if SGData['SGLaue'] in ['-1',]:
612	return [pfx+'A0',pfx+'A1',pfx+'A2',pfx+'A3',pfx+'A4',pfx+'A5']
613	elif SGData['SGLaue'] in ['2/m',]:
614	if SGData['SGUniq'] == 'a':
615	return [pfx+'A0',pfx+'A1',pfx+'A2',pfx+'A3']
616	elif SGData['SGUniq'] == 'b':
617	return [pfx+'A0',pfx+'A1',pfx+'A2',pfx+'A4']
618	else:
619	return [pfx+'A0',pfx+'A1',pfx+'A2',pfx+'A5']
620	elif SGData['SGLaue'] in ['mmm',]:
621	return [pfx+'A0',pfx+'A1',pfx+'A2']
622	elif SGData['SGLaue'] in ['4/m','4/mmm']:
623	return [pfx+'A0',pfx+'A2']
624	elif SGData['SGLaue'] in ['6/m','6/mmm','3m1', '31m', '3']:
625	return [pfx+'A0',pfx+'A2']
626	elif SGData['SGLaue'] in ['3R', '3mR']:
627	return [pfx+'A0',pfx+'A3']
628	elif SGData['SGLaue'] in ['m3m','m3']:
629	return [pfx+'A0',]
630
631	################################################################################
632	##### Rigid Body Models and not General.get('doPawley')
633	################################################################################
634
635	def GetRigidBodyModels(rigidbodyDict,Print=True,pFile=None):
636	'needs a doc string'
637
638	def PrintResRBModel(RBModel):
639	atNames = RBModel['atNames']
640	rbRef = RBModel['rbRef']
641	rbSeq = RBModel['rbSeq']
642	print >>pFile,'Residue RB name: ',RBModel['RBname'],' no.atoms: ',len(RBModel['rbTypes']), \
643	'No. times used: ',RBModel['useCount']
644	print >>pFile,' At name x y z'
645	for name,xyz in zip(atNames,RBModel['rbXYZ']):
646	print >>pFile,' %8s %10.4f %10.4f %10.4f'%(name,xyz[0],xyz[1],xyz[2])
647	print >>pFile,'Orientation defined by:',atNames[rbRef[0]],' -> ',atNames[rbRef[1]], \
648	' & ',atNames[rbRef[0]],' -> ',atNames[rbRef[2]]
649	if rbSeq:
650	for i,rbseq in enumerate(rbSeq):
651	print >>pFile,'Torsion sequence ',i,' Bond: '+atNames[rbseq[0]],' - ', \
652	atNames[rbseq[1]],' riding: ',[atNames[i] for i in rbseq[3]]
653
654	def PrintVecRBModel(RBModel):
655	rbRef = RBModel['rbRef']
656	atTypes = RBModel['rbTypes']
657	print >>pFile,'Vector RB name: ',RBModel['RBname'],' no.atoms: ',len(RBModel['rbTypes']), \
658	'No. times used: ',RBModel['useCount']
659	for i in range(len(RBModel['VectMag'])):
660	print >>pFile,'Vector no.: ',i,' Magnitude: ', \
661	'%8.4f'%(RBModel['VectMag'][i]),' Refine? ',RBModel['VectRef'][i]
662	print >>pFile,' No. Type vx vy vz'
663	for j,[name,xyz] in enumerate(zip(atTypes,RBModel['rbVect'][i])):
664	print >>pFile,' %d %2s %10.4f %10.4f %10.4f'%(j,name,xyz[0],xyz[1],xyz[2])
665	print >>pFile,' No. Type x y z'
666	for i,[name,xyz] in enumerate(zip(atTypes,RBModel['rbXYZ'])):
667	print >>pFile,' %d %2s %10.4f %10.4f %10.4f'%(i,name,xyz[0],xyz[1],xyz[2])
668	print >>pFile,'Orientation defined by: atom ',rbRef[0],' -> atom ',rbRef[1], \
669	' & atom ',rbRef[0],' -> atom ',rbRef[2]
670	rbVary = []
671	rbDict = {}
672	rbIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[]})
673	if len(rbIds['Vector']):
674	for irb,item in enumerate(rbIds['Vector']):
675	if rigidbodyDict['Vector'][item]['useCount']:
676	RBmags = rigidbodyDict['Vector'][item]['VectMag']
677	RBrefs = rigidbodyDict['Vector'][item]['VectRef']
678	for i,[mag,ref] in enumerate(zip(RBmags,RBrefs)):
679	pid = '::RBV;'+str(i)+':'+str(irb)
680	rbDict[pid] = mag
681	if ref:
682	rbVary.append(pid)
683	if Print:
684	print >>pFile,'\nVector rigid body model:'
685	PrintVecRBModel(rigidbodyDict['Vector'][item])
686	if len(rbIds['Residue']):
687	for item in rbIds['Residue']:
688	if rigidbodyDict['Residue'][item]['useCount']:
689	if Print:
690	print >>pFile,'\nResidue rigid body model:'
691	PrintResRBModel(rigidbodyDict['Residue'][item])
692	return rbVary,rbDict
693
694	def SetRigidBodyModels(parmDict,sigDict,rigidbodyDict,pFile=None):
695	'needs a doc string'
696
697	def PrintRBVectandSig(VectRB,VectSig):
698	print >>pFile,'\n Rigid body vector magnitudes for '+VectRB['RBname']+':'
699	namstr = ' names :'
700	valstr = ' values:'
701	sigstr = ' esds :'
702	for i,[val,sig] in enumerate(zip(VectRB['VectMag'],VectSig)):
703	namstr += '%12s'%('Vect '+str(i))
704	valstr += '%12.4f'%(val)
705	if sig:
706	sigstr += '%12.4f'%(sig)
707	else:
708	sigstr += 12*' '
709	print >>pFile,namstr
710	print >>pFile,valstr
711	print >>pFile,sigstr
712
713	RBIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[]}) #these are lists of rbIds
714	if not RBIds['Vector']:
715	return
716	for irb,item in enumerate(RBIds['Vector']):
717	if rigidbodyDict['Vector'][item]['useCount']:
718	VectSig = []
719	RBmags = rigidbodyDict['Vector'][item]['VectMag']
720	for i,mag in enumerate(RBmags):
721	name = '::RBV;'+str(i)+':'+str(irb)
722	mag = parmDict[name]
723	if name in sigDict:
724	VectSig.append(sigDict[name])
725	PrintRBVectandSig(rigidbodyDict['Vector'][item],VectSig)
726
727	################################################################################
728	##### Phase data
729	################################################################################
730
731	def GetPhaseData(PhaseData,RestraintDict={},rbIds={},Print=True,pFile=None):
732	'needs a doc string'
733
734	def PrintFFtable(FFtable):
735	print >>pFile,'\n X-ray scattering factors:'
736	print >>pFile,' Symbol fa fb fc'
737	print >>pFile,99*'-'
738	for Ename in FFtable:
739	ffdata = FFtable[Ename]
740	fa = ffdata['fa']
741	fb = ffdata['fb']
742	print >>pFile,' %8s %9.5f %9.5f %9.5f %9.5f %9.5f %9.5f %9.5f %9.5f %9.5f' % \
743	(Ename.ljust(8),fa[0],fa[1],fa[2],fa[3],fb[0],fb[1],fb[2],fb[3],ffdata['fc'])
744
745	def PrintBLtable(BLtable):
746	print >>pFile,'\n Neutron scattering factors:'
747	print >>pFile,' Symbol isotope mass b resonant terms'
748	print >>pFile,99*'-'
749	for Ename in BLtable:
750	bldata = BLtable[Ename]
751	isotope = bldata[0]
752	mass = bldata[1]['Mass']
753	if 'SL' in bldata[1]:
754	blen = bldata[1]['SL'][0]
755	bres = []
756	else:
757	blen = 0
758	bres = bldata[1]['BW-LS']
759	line = ' %8s%11s %10.3f %8.3f'%(Ename.ljust(8),isotope.center(11),mass,blen)
760	for item in bres:
761	line += '%10.5g'%(item)
762	print >>pFile,line
763
764	def PrintRBObjects(resRBData,vecRBData):
765
766	def PrintRBThermals():
767	tlstr = ['11','22','33','12','13','23']
768	sstr = ['12','13','21','23','31','32','AA','BB']
769	TLS = RB['ThermalMotion'][1]
770	TLSvar = RB['ThermalMotion'][2]
771	if 'T' in RB['ThermalMotion'][0]:
772	print >>pFile,'TLS data'
773	text = ''
774	for i in range(6):
775	text += 'T'+tlstr[i]+' %8.4f %s '%(TLS[i],str(TLSvar[i])[0])
776	print >>pFile,text
777	if 'L' in RB['ThermalMotion'][0]:
778	text = ''
779	for i in range(6,12):
780	text += 'L'+tlstr[i-6]+' %8.2f %s '%(TLS[i],str(TLSvar[i])[0])
781	print >>pFile,text
782	if 'S' in RB['ThermalMotion'][0]:
783	text = ''
784	for i in range(12,20):
785	text += 'S'+sstr[i-12]+' %8.3f %s '%(TLS[i],str(TLSvar[i])[0])
786	print >>pFile,text
787	if 'U' in RB['ThermalMotion'][0]:
788	print >>pFile,'Uiso data'
789	text = 'Uiso'+' %10.3f %s'%(TLS[0],str(TLSvar[0])[0])
790
791	if len(resRBData):
792	for RB in resRBData:
793	Oxyz = RB['Orig'][0]
794	Qrijk = RB['Orient'][0]
795	Angle = 2.0*acosd(Qrijk[0])
796	print >>pFile,'\nRBObject ',RB['RBname'],' at ', \
797	'%10.4f %10.4f %10.4f'%(Oxyz[0],Oxyz[1],Oxyz[2]),' Refine?',RB['Orig'][1]
798	print >>pFile,'Orientation angle,vector:', \
799	'%10.3f %10.4f %10.4f %10.4f'%(Angle,Qrijk[1],Qrijk[2],Qrijk[3]),' Refine? ',RB['Orient'][1]
800	Torsions = RB['Torsions']
801	if len(Torsions):
802	text = 'Torsions: '
803	for torsion in Torsions:
804	text += '%10.4f Refine? %s'%(torsion[0],torsion[1])
805	print >>pFile,text
806	PrintRBThermals()
807	if len(vecRBData):
808	for RB in vecRBData:
809	Oxyz = RB['Orig'][0]
810	Qrijk = RB['Orient'][0]
811	Angle = 2.0*acosd(Qrijk[0])
812	print >>pFile,'\nRBObject ',RB['RBname'],' at ', \
813	'%10.4f %10.4f %10.4f'%(Oxyz[0],Oxyz[1],Oxyz[2]),' Refine?',RB['Orig'][1]
814	print >>pFile,'Orientation angle,vector:', \
815	'%10.3f %10.4f %10.4f %10.4f'%(Angle,Qrijk[1],Qrijk[2],Qrijk[3]),' Refine? ',RB['Orient'][1]
816	PrintRBThermals()
817
818	def PrintAtoms(General,Atoms):
819	cx,ct,cs,cia = General['AtomPtrs']
820	print >>pFile,'\n Atoms:'
821	line = ' name type refine? x y z '+ \
822	' frac site sym mult I/A Uiso U11 U22 U33 U12 U13 U23'
823	if General['Type'] == 'magnetic':
824	line += ' Mx My Mz'
825	elif General['Type'] == 'macromolecular':
826	line = ' res no residue chain'+line
827	print >>pFile,line
828	if General['Type'] == 'nuclear':
829	print >>pFile,135*'-'
830	for i,at in enumerate(Atoms):
831	line = '%7s'%(at[ct-1])+'%7s'%(at[ct])+'%7s'%(at[ct+1])+'%10.5f'%(at[cx])+'%10.5f'%(at[cx+1])+ \
832	'%10.5f'%(at[cx+2])+'%8.3f'%(at[cx+3])+'%7s'%(at[cs])+'%5d'%(at[cs+1])+'%5s'%(at[cia])
833	if at[cia] == 'I':
834	line += '%8.5f'%(at[cia+1])+48*' '
835	else:
836	line += 8*' '
837	for j in range(6):
838	line += '%8.5f'%(at[cia+2+j])
839	print >>pFile,line
840	elif General['Type'] == 'macromolecular':
841	print >>pFile,135*'-'
842	for i,at in enumerate(Atoms):
843	line = '%7s'%(at[0])+'%7s'%(at[1])+'%7s'%(at[2])+'%7s'%(at[ct-1])+'%7s'%(at[ct])+'%7s'%(at[ct+1])+'%10.5f'%(at[cx])+'%10.5f'%(at[cx+1])+ \
844	'%10.5f'%(at[cx+2])+'%8.3f'%(at[cx+3])+'%7s'%(at[cs])+'%5d'%(at[cs+1])+'%5s'%(at[cia])
845	if at[cia] == 'I':
846	line += '%8.4f'%(at[cia+1])+48*' '
847	else:
848	line += 8*' '
849	for j in range(6):
850	line += '%8.4f'%(at[cia+2+j])
851	print >>pFile,line
852
853	def PrintTexture(textureData):
854	topstr = '\n Spherical harmonics texture: Order:' + \
855	str(textureData['Order'])
856	if textureData['Order']:
857	print >>pFile,topstr+' Refine? '+str(textureData['SH Coeff'][0])
858	else:
859	print >>pFile,topstr
860	return
861	names = ['omega','chi','phi']
862	line = '\n'
863	for name in names:
864	line += ' SH '+name+':'+'%12.4f'%(textureData['Sample '+name][1])+' Refine? '+str(textureData['Sample '+name][0])
865	print >>pFile,line
866	print >>pFile,'\n Texture coefficients:'
867	ptlbls = ' names :'
868	ptstr = ' values:'
869	SHcoeff = textureData['SH Coeff'][1]
870	for item in SHcoeff:
871	ptlbls += '%12s'%(item)
872	ptstr += '%12.4f'%(SHcoeff[item])
873	print >>pFile,ptlbls
874	print >>pFile,ptstr
875
876	def MakeRBParms(rbKey,phaseVary,phaseDict):
877	rbid = str(rbids.index(RB['RBId']))
878	pfxRB = pfx+'RB'+rbKey+'P'
879	pstr = ['x','y','z']
880	ostr = ['a','i','j','k']
881	for i in range(3):
882	name = pfxRB+pstr[i]+':'+str(iRB)+':'+rbid
883	phaseDict[name] = RB['Orig'][0][i]
884	if RB['Orig'][1]:
885	phaseVary += [name,]
886	pfxRB = pfx+'RB'+rbKey+'O'
887	for i in range(4):
888	name = pfxRB+ostr[i]+':'+str(iRB)+':'+rbid
889	phaseDict[name] = RB['Orient'][0][i]
890	if RB['Orient'][1] == 'AV' and i:
891	phaseVary += [name,]
892	elif RB['Orient'][1] == 'A' and not i:
893	phaseVary += [name,]
894
895	def MakeRBThermals(rbKey,phaseVary,phaseDict):
896	rbid = str(rbids.index(RB['RBId']))
897	tlstr = ['11','22','33','12','13','23']
898	sstr = ['12','13','21','23','31','32','AA','BB']
899	if 'T' in RB['ThermalMotion'][0]:
900	pfxRB = pfx+'RB'+rbKey+'T'
901	for i in range(6):
902	name = pfxRB+tlstr[i]+':'+str(iRB)+':'+rbid
903	phaseDict[name] = RB['ThermalMotion'][1][i]
904	if RB['ThermalMotion'][2][i]:
905	phaseVary += [name,]
906	if 'L' in RB['ThermalMotion'][0]:
907	pfxRB = pfx+'RB'+rbKey+'L'
908	for i in range(6):
909	name = pfxRB+tlstr[i]+':'+str(iRB)+':'+rbid
910	phaseDict[name] = RB['ThermalMotion'][1][i+6]
911	if RB['ThermalMotion'][2][i+6]:
912	phaseVary += [name,]
913	if 'S' in RB['ThermalMotion'][0]:
914	pfxRB = pfx+'RB'+rbKey+'S'
915	for i in range(8):
916	name = pfxRB+sstr[i]+':'+str(iRB)+':'+rbid
917	phaseDict[name] = RB['ThermalMotion'][1][i+12]
918	if RB['ThermalMotion'][2][i+12]:
919	phaseVary += [name,]
920	if 'U' in RB['ThermalMotion'][0]:
921	name = pfx+'RB'+rbKey+'U:'+str(iRB)+':'+rbid
922	phaseDict[name] = RB['ThermalMotion'][1][0]
923	if RB['ThermalMotion'][2][0]:
924	phaseVary += [name,]
925
926	def MakeRBTorsions(rbKey,phaseVary,phaseDict):
927	rbid = str(rbids.index(RB['RBId']))
928	pfxRB = pfx+'RB'+rbKey+'Tr;'
929	for i,torsion in enumerate(RB['Torsions']):
930	name = pfxRB+str(i)+':'+str(iRB)+':'+rbid
931	phaseDict[name] = torsion[0]
932	if torsion[1]:
933	phaseVary += [name,]
934
935	if Print:
936	print >>pFile,'\n Phases:'
937	phaseVary = []
938	phaseDict = {}
939	phaseConstr = {}
940	pawleyLookup = {}
941	FFtables = {} #scattering factors - xrays
942	BLtables = {} # neutrons
943	Natoms = {}
944	AtMults = {}
945	AtIA = {}
946	shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
947	SamSym = dict(zip(shModels,['0','-1','2/m','mmm']))
948	atomIndx = {}
949	for name in PhaseData:
950	General = PhaseData[name]['General']
951	pId = PhaseData[name]['pId']
952	pfx = str(pId)+'::'
953	FFtable = G2el.GetFFtable(General['AtomTypes'])
954	BLtable = G2el.GetBLtable(General)
955	FFtables.update(FFtable)
956	BLtables.update(BLtable)
957	Atoms = PhaseData[name]['Atoms']
958	AtLookup = G2mth.FillAtomLookUp(Atoms)
959	PawleyRef = PhaseData[name].get('Pawley ref',[])
960	SGData = General['SGData']
961	SGtext,SGtable = G2spc.SGPrint(SGData)
962	cell = General['Cell']
963	A = G2lat.cell2A(cell[1:7])
964	phaseDict.update({pfx+'A0':A[0],pfx+'A1':A[1],pfx+'A2':A[2],
965	pfx+'A3':A[3],pfx+'A4':A[4],pfx+'A5':A[5],pfx+'Vol':G2lat.calc_V(A)})
966	if cell[0]:
967	phaseVary += cellVary(pfx,SGData)
968	resRBData = PhaseData[name]['RBModels'].get('Residue',[])
969	if resRBData:
970	rbids = rbIds['Residue'] #NB: used in the MakeRB routines
971	for iRB,RB in enumerate(resRBData):
972	MakeRBParms('R',phaseVary,phaseDict)
973	MakeRBThermals('R',phaseVary,phaseDict)
974	MakeRBTorsions('R',phaseVary,phaseDict)
975
976	vecRBData = PhaseData[name]['RBModels'].get('Vector',[])
977	if vecRBData:
978	rbids = rbIds['Vector'] #NB: used in the MakeRB routines
979	for iRB,RB in enumerate(vecRBData):
980	MakeRBParms('V',phaseVary,phaseDict)
981	MakeRBThermals('V',phaseVary,phaseDict)
982
983	Natoms[pfx] = 0
984	if Atoms and not General.get('doPawley'):
985	cx,ct,cs,cia = General['AtomPtrs']
986	if General['Type'] in ['nuclear','macromolecular']:
987	Natoms[pfx] = len(Atoms)
988	for i,at in enumerate(Atoms):
989	atomIndx[at[-1]] = [pfx,i] #lookup table for restraints
990	phaseDict.update({pfx+'Atype:'+str(i):at[ct],pfx+'Afrac:'+str(i):at[cx+3],pfx+'Amul:'+str(i):at[cs+1],
991	pfx+'Ax:'+str(i):at[cx],pfx+'Ay:'+str(i):at[cx+1],pfx+'Az:'+str(i):at[cx+2],
992	pfx+'dAx:'+str(i):0.,pfx+'dAy:'+str(i):0.,pfx+'dAz:'+str(i):0., #refined shifts for x,y,z
993	pfx+'AI/A:'+str(i):at[cia],})
994	if at[cia] == 'I':
995	phaseDict[pfx+'AUiso:'+str(i)] = at[cia+1]
996	else:
997	phaseDict.update({pfx+'AU11:'+str(i):at[cia+2],pfx+'AU22:'+str(i):at[cia+3],pfx+'AU33:'+str(i):at[cia+4],
998	pfx+'AU12:'+str(i):at[cia+5],pfx+'AU13:'+str(i):at[cia+6],pfx+'AU23:'+str(i):at[cia+7]})
999	if 'F' in at[ct+1]:
1000	phaseVary.append(pfx+'Afrac:'+str(i))
1001	if 'X' in at[ct+1]:
1002	try: #patch for sytsym name changes
1003	xId,xCoef = G2spc.GetCSxinel(at[cs])
1004	except KeyError:
1005	Sytsym = G2spc.SytSym(at[cx:cx+3],SGData)[0]
1006	at[cs] = Sytsym
1007	xId,xCoef = G2spc.GetCSxinel(at[cs])
1008	xId,xCoef = G2spc.GetCSxinel(at[cs])
1009	names = [pfx+'dAx:'+str(i),pfx+'dAy:'+str(i),pfx+'dAz:'+str(i)]
1010	equivs = [[],[],[]]
1011	for j in range(3):
1012	if xId[j] > 0:
1013	phaseVary.append(names[j])
1014	equivs[xId[j]-1].append([names[j],xCoef[j]])
1015	for equiv in equivs:
1016	if len(equiv) > 1:
1017	name = equiv[0][0]
1018	coef = equiv[0][1]
1019	for eqv in equiv[1:]:
1020	eqv[1] /= coef
1021	G2mv.StoreEquivalence(name,(eqv,))
1022	if 'U' in at[ct+1]:
1023	if at[cia] == 'I':
1024	phaseVary.append(pfx+'AUiso:'+str(i))
1025	else:
1026	try: #patch for sytsym name changes
1027	uId,uCoef = G2spc.GetCSuinel(at[cs])[:2]
1028	except KeyError:
1029	Sytsym = G2spc.SytSym(at[cx:cx+3],SGData)[0]
1030	at[cs] = Sytsym
1031	uId,uCoef = G2spc.GetCSuinel(at[cs])[:2]
1032	uId,uCoef = G2spc.GetCSuinel(at[cs])[:2]
1033	names = [pfx+'AU11:'+str(i),pfx+'AU22:'+str(i),pfx+'AU33:'+str(i),
1034	pfx+'AU12:'+str(i),pfx+'AU13:'+str(i),pfx+'AU23:'+str(i)]
1035	equivs = [[],[],[],[],[],[]]
1036	for j in range(6):
1037	if uId[j] > 0:
1038	phaseVary.append(names[j])
1039	equivs[uId[j]-1].append([names[j],uCoef[j]])
1040	for equiv in equivs:
1041	if len(equiv) > 1:
1042	name = equiv[0][0]
1043	coef = equiv[0][1]
1044	for eqv in equiv[1:]:
1045	eqv[1] /= coef
1046	G2mv.StoreEquivalence(name,equiv[1:])
1047	# elif General['Type'] == 'magnetic':
1048	# elif General['Type'] == 'macromolecular':
1049	textureData = General['SH Texture']
1050	if textureData['Order']:
1051	phaseDict[pfx+'SHorder'] = textureData['Order']
1052	phaseDict[pfx+'SHmodel'] = SamSym[textureData['Model']]
1053	for item in ['omega','chi','phi']:
1054	phaseDict[pfx+'SH '+item] = textureData['Sample '+item][1]
1055	if textureData['Sample '+item][0]:
1056	phaseVary.append(pfx+'SH '+item)
1057	for item in textureData['SH Coeff'][1]:
1058	phaseDict[pfx+item] = textureData['SH Coeff'][1][item]
1059	if textureData['SH Coeff'][0]:
1060	phaseVary.append(pfx+item)
1061
1062	if Print:
1063	print >>pFile,'\n Phase name: ',General['Name']
1064	print >>pFile,135*'-'
1065	PrintFFtable(FFtable)
1066	PrintBLtable(BLtable)
1067	print >>pFile,''
1068	for line in SGtext: print >>pFile,line
1069	if len(SGtable):
1070	for item in SGtable:
1071	line = ' %s '%(item)
1072	print >>pFile,line
1073	else:
1074	print >>pFile,' ( 1) %s'%(SGtable[0])
1075	PrintRBObjects(resRBData,vecRBData)
1076	PrintAtoms(General,Atoms)
1077	print >>pFile,'\n Unit cell: a =','%.5f'%(cell[1]),' b =','%.5f'%(cell[2]),' c =','%.5f'%(cell[3]), \
1078	' alpha =','%.3f'%(cell[4]),' beta =','%.3f'%(cell[5]),' gamma =', \
1079	'%.3f'%(cell[6]),' volume =','%.3f'%(cell[7]),' Refine?',cell[0]
1080	PrintTexture(textureData)
1081	if name in RestraintDict:
1082	PrintRestraints(cell[1:7],SGData,General['AtomPtrs'],Atoms,AtLookup,
1083	textureData,RestraintDict[name],pFile)
1084
1085	elif PawleyRef:
1086	pawleyVary = []
1087	for i,refl in enumerate(PawleyRef):
1088	phaseDict[pfx+'PWLref:'+str(i)] = refl[6]
1089	pawleyLookup[pfx+'%d,%d,%d'%(refl[0],refl[1],refl[2])] = i
1090	if refl[5]:
1091	pawleyVary.append(pfx+'PWLref:'+str(i))
1092	GetPawleyConstr(SGData['SGLaue'],PawleyRef,pawleyVary) #does G2mv.StoreEquivalence
1093	phaseVary += pawleyVary
1094
1095	return Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables
1096
1097	def cellFill(pfx,SGData,parmDict,sigDict):
1098	'''Returns the filled-out reciprocal cell (A) terms and their uncertainties
1099	from the parameter and sig dictionaries.
1100
1101	:param str pfx: parameter prefix ("n::", where n is a phase number)
1102	:param dict SGdata: a symmetry object
1103	:param dict parmDict: a dictionary of parameters
1104	:param dict sigDict: a dictionary of uncertainties on parameters
1105
1106	:returns: A,sigA where each is a list of six terms with the A terms
1107	'''
1108	if SGData['SGLaue'] in ['-1',]:
1109	A = [parmDict[pfx+'A0'],parmDict[pfx+'A1'],parmDict[pfx+'A2'],
1110	parmDict[pfx+'A3'],parmDict[pfx+'A4'],parmDict[pfx+'A5']]
1111	elif SGData['SGLaue'] in ['2/m',]:
1112	if SGData['SGUniq'] == 'a':
1113	A = [parmDict[pfx+'A0'],parmDict[pfx+'A1'],parmDict[pfx+'A2'],
1114	parmDict[pfx+'A3'],0,0]
1115	elif SGData['SGUniq'] == 'b':
1116	A = [parmDict[pfx+'A0'],parmDict[pfx+'A1'],parmDict[pfx+'A2'],
1117	0,parmDict[pfx+'A4'],0]
1118	else:
1119	A = [parmDict[pfx+'A0'],parmDict[pfx+'A1'],parmDict[pfx+'A2'],
1120	0,0,parmDict[pfx+'A5']]
1121	elif SGData['SGLaue'] in ['mmm',]:
1122	A = [parmDict[pfx+'A0'],parmDict[pfx+'A1'],parmDict[pfx+'A2'],0,0,0]
1123	elif SGData['SGLaue'] in ['4/m','4/mmm']:
1124	A = [parmDict[pfx+'A0'],parmDict[pfx+'A0'],parmDict[pfx+'A2'],0,0,0]
1125	elif SGData['SGLaue'] in ['6/m','6/mmm','3m1', '31m', '3']:
1126	A = [parmDict[pfx+'A0'],parmDict[pfx+'A0'],parmDict[pfx+'A2'],
1127	parmDict[pfx+'A0'],0,0]
1128	elif SGData['SGLaue'] in ['3R', '3mR']:
1129	A = [parmDict[pfx+'A0'],parmDict[pfx+'A0'],parmDict[pfx+'A0'],
1130	parmDict[pfx+'A3'],parmDict[pfx+'A3'],parmDict[pfx+'A3']]
1131	elif SGData['SGLaue'] in ['m3m','m3']:
1132	A = [parmDict[pfx+'A0'],parmDict[pfx+'A0'],parmDict[pfx+'A0'],0,0,0]
1133
1134	try:
1135	if SGData['SGLaue'] in ['-1',]:
1136	sigA = [sigDict[pfx+'A0'],sigDict[pfx+'A1'],sigDict[pfx+'A2'],
1137	sigDict[pfx+'A3'],sigDict[pfx+'A4'],sigDict[pfx+'A5']]
1138	elif SGData['SGLaue'] in ['2/m',]:
1139	if SGData['SGUniq'] == 'a':
1140	sigA = [sigDict[pfx+'A0'],sigDict[pfx+'A1'],sigDict[pfx+'A2'],
1141	sigDict[pfx+'A3'],0,0]
1142	elif SGData['SGUniq'] == 'b':
1143	sigA = [sigDict[pfx+'A0'],sigDict[pfx+'A1'],sigDict[pfx+'A2'],
1144	0,sigDict[pfx+'A4'],0]
1145	else:
1146	sigA = [sigDict[pfx+'A0'],sigDict[pfx+'A1'],sigDict[pfx+'A2'],
1147	0,0,sigDict[pfx+'A5']]
1148	elif SGData['SGLaue'] in ['mmm',]:
1149	sigA = [sigDict[pfx+'A0'],sigDict[pfx+'A1'],sigDict[pfx+'A2'],0,0,0]
1150	elif SGData['SGLaue'] in ['4/m','4/mmm']:
1151	sigA = [sigDict[pfx+'A0'],0,sigDict[pfx+'A2'],0,0,0]
1152	elif SGData['SGLaue'] in ['6/m','6/mmm','3m1', '31m', '3']:
1153	sigA = [sigDict[pfx+'A0'],0,sigDict[pfx+'A2'],0,0,0]
1154	elif SGData['SGLaue'] in ['3R', '3mR']:
1155	sigA = [sigDict[pfx+'A0'],0,0,sigDict[pfx+'A3'],0,0]
1156	elif SGData['SGLaue'] in ['m3m','m3']:
1157	sigA = [sigDict[pfx+'A0'],0,0,0,0,0]
1158	except KeyError:
1159	sigA = [0,0,0,0,0,0]
1160
1161	return A,sigA
1162
1163	def PrintRestraints(cell,SGData,AtPtrs,Atoms,AtLookup,textureData,phaseRest,pFile):
1164	'needs a doc string'
1165	if phaseRest:
1166	Amat = G2lat.cell2AB(cell)[0]
1167	cx,ct,cs = AtPtrs[:3]
1168	names = [['Bond','Bonds'],['Angle','Angles'],['Plane','Planes'],
1169	['Chiral','Volumes'],['Torsion','Torsions'],['Rama','Ramas'],
1170	['ChemComp','Sites'],['Texture','HKLs']]
1171	for name,rest in names:
1172	itemRest = phaseRest[name]
1173	if itemRest[rest] and itemRest['Use']:
1174	print >>pFile,'\n %s %10.3f Use: %s'%(name+' restraint weight factor',itemRest['wtFactor'],str(itemRest['Use']))
1175	if name in ['Bond','Angle','Plane','Chiral']:
1176	print >>pFile,' calc obs sig delt/sig atoms(symOp): '
1177	for indx,ops,obs,esd in itemRest[rest]:
1178	try:
1179	AtNames = G2mth.GetAtomItemsById(Atoms,AtLookup,indx,ct-1)
1180	AtName = ''
1181	for i,Aname in enumerate(AtNames):
1182	AtName += Aname
1183	if ops[i] == '1':
1184	AtName += '-'
1185	else:
1186	AtName += '+('+ops[i]+')-'
1187	XYZ = np.array(G2mth.GetAtomItemsById(Atoms,AtLookup,indx,cx,3))
1188	XYZ = G2mth.getSyXYZ(XYZ,ops,SGData)
1189	if name == 'Bond':
1190	calc = G2mth.getRestDist(XYZ,Amat)
1191	elif name == 'Angle':
1192	calc = G2mth.getRestAngle(XYZ,Amat)
1193	elif name == 'Plane':
1194	calc = G2mth.getRestPlane(XYZ,Amat)
1195	elif name == 'Chiral':
1196	calc = G2mth.getRestChiral(XYZ,Amat)
1197	print >>pFile,' %9.3f %9.3f %8.3f %8.3f %s'%(calc,obs,esd,(obs-calc)/esd,AtName[:-1])
1198	except KeyError:
1199	del itemRest[rest]
1200	elif name in ['Torsion','Rama']:
1201	print >>pFile,' atoms(symOp) calc obs sig delt/sig torsions: '
1202	coeffDict = itemRest['Coeff']
1203	for indx,ops,cofName,esd in itemRest[rest]:
1204	AtNames = G2mth.GetAtomItemsById(Atoms,AtLookup,indx,ct-1)
1205	AtName = ''
1206	for i,Aname in enumerate(AtNames):
1207	AtName += Aname+'+('+ops[i]+')-'
1208	XYZ = np.array(G2mth.GetAtomItemsById(Atoms,AtLookup,indx,cx,3))
1209	XYZ = G2mth.getSyXYZ(XYZ,ops,SGData)
1210	if name == 'Torsion':
1211	tor = G2mth.getRestTorsion(XYZ,Amat)
1212	restr,calc = G2mth.calcTorsionEnergy(tor,coeffDict[cofName])
1213	print >>pFile,' %8.3f %8.3f %.3f %8.3f %8.3f %s'%(calc,obs,esd,(obs-calc)/esd,tor,AtName[:-1])
1214	else:
1215	phi,psi = G2mth.getRestRama(XYZ,Amat)
1216	restr,calc = G2mth.calcRamaEnergy(phi,psi,coeffDict[cofName])
1217	print >>pFile,' %8.3f %8.3f %8.3f %8.3f %8.3f %8.3f %s'%(calc,obs,esd,(obs-calc)/esd,phi,psi,AtName[:-1])
1218	elif name == 'ChemComp':
1219	print >>pFile,' atoms mul*frac factor prod'
1220	for indx,factors,obs,esd in itemRest[rest]:
1221	try:
1222	atoms = G2mth.GetAtomItemsById(Atoms,AtLookup,indx,ct-1)
1223	mul = np.array(G2mth.GetAtomItemsById(Atoms,AtLookup,indx,cs+1))
1224	frac = np.array(G2mth.GetAtomItemsById(Atoms,AtLookup,indx,cs-1))
1225	mulfrac = mul*frac
1226	calcs = mulfracfactors
1227	for iatm,[atom,mf,fr,clc] in enumerate(zip(atoms,mulfrac,factors,calcs)):
1228	print >>pFile,' %10s %8.3f %8.3f %8.3f'%(atom,mf,fr,clc)
1229	print >>pFile,' Sum: calc: %8.3f obs: %8.3f esd: %8.3f'%(np.sum(calcs),obs,esd)
1230	except KeyError:
1231	del itemRest[rest]
1232	elif name == 'Texture' and textureData['Order']:
1233	Start = False
1234	SHCoef = textureData['SH Coeff'][1]
1235	shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
1236	SamSym = dict(zip(shModels,['0','-1','2/m','mmm']))
1237	print ' HKL grid neg esd sum neg num neg use unit? unit esd '
1238	for hkl,grid,esd1,ifesd2,esd2 in itemRest[rest]:
1239	PH = np.array(hkl)
1240	phi,beta = G2lat.CrsAng(np.array(hkl),cell,SGData)
1241	ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1242	R,P,Z = G2mth.getRestPolefig(ODFln,SamSym[textureData['Model']],grid)
1243	Z = ma.masked_greater(Z,0.0)
1244	num = ma.count(Z)
1245	sum = 0
1246	if num:
1247	sum = np.sum(Z)
1248	print ' %d %d %d %d %8.3f %8.3f %8d %s %8.3f'%(hkl[0],hkl[1],hkl[2],grid,esd1,sum,num,str(ifesd2),esd2)
1249
1250	def getCellEsd(pfx,SGData,A,covData):
1251	'needs a doc string'
1252	dpr = 180./np.pi
1253	rVsq = G2lat.calc_rVsq(A)
1254	G,g = G2lat.A2Gmat(A) #get recip. & real metric tensors
1255	cell = np.array(G2lat.Gmat2cell(g)) #real cell
1256	cellst = np.array(G2lat.Gmat2cell(G)) #recip. cell
1257	scos = cosd(cellst[3:6])
1258	ssin = sind(cellst[3:6])
1259	scot = scos/ssin
1260	rcos = cosd(cell[3:6])
1261	rsin = sind(cell[3:6])
1262	rcot = rcos/rsin
1263	RMnames = [pfx+'A0',pfx+'A1',pfx+'A2',pfx+'A3',pfx+'A4',pfx+'A5']
1264	varyList = covData['varyList']
1265	covMatrix = covData['covMatrix']
1266	vcov = G2mth.getVCov(RMnames,varyList,covMatrix)
1267	Ax = np.array(A)
1268	Ax[3:] /= 2.
1269	drVdA = np.array([Ax[1]Ax[2]-Ax[5]2,Ax[0]Ax[2]-Ax[4]*2,Ax[0]Ax[1]-Ax[3]**2,
1270	Ax[4]Ax[5]-Ax[2]Ax[3],Ax[3]Ax[5]-Ax[1]Ax[4],Ax[3]Ax[4]-Ax[0]Ax[5]])
1271	srcvlsq = np.inner(drVdA,np.inner(vcov,drVdA.T))
1272	Vol = 1/np.sqrt(rVsq)
1273	sigVol = Vol*3np.sqrt(srcvlsq)/2.
1274	R123 = Ax[0]Ax[1]Ax[2]
1275	dsasdg = np.zeros((3,6))
1276	dadg = np.zeros((6,6))
1277	for i0 in range(3): #0 1 2
1278	i1 = (i0+1)%3 #1 2 0
1279	i2 = (i1+1)%3 #2 0 1
1280	i3 = 5-i2 #3 5 4
1281	i4 = 5-i1 #4 3 5
1282	i5 = 5-i0 #5 4 3
1283	dsasdg[i0][i1] = 0.5scot[i0]scos[i0]/Ax[i1]
1284	dsasdg[i0][i2] = 0.5scot[i0]scos[i0]/Ax[i2]
1285	dsasdg[i0][i5] = -scot[i0]/np.sqrt(Ax[i1]*Ax[i2])
1286	denmsq = Ax[i0](R123-Ax[i1]Ax[i4]*2-Ax[i2]Ax[i3]*2+(Ax[i4]Ax[i3])**2)
1287	denom = np.sqrt(denmsq)
1288	dadg[i5][i0] = -Ax[i5]/denom-rcos[i0]/denmsq(R123-0.5Ax[i1]Ax[i4]2-0.5Ax[i2]Ax[i3]*2)
1289	dadg[i5][i1] = -0.5rcos[i0]/denmsq(Ax[i0]*2Ax[i2]-Ax[i0]Ax[i4]*2)
1290	dadg[i5][i2] = -0.5rcos[i0]/denmsq(Ax[i0]*2Ax[i1]-Ax[i0]Ax[i3]*2)
1291	dadg[i5][i3] = Ax[i4]/denom+rcos[i0]/denmsq(Ax[i0]Ax[i2]Ax[i3]-Ax[i3]Ax[i4]**2)
1292	dadg[i5][i4] = Ax[i3]/denom+rcos[i0]/denmsq(Ax[i0]Ax[i1]Ax[i4]-Ax[i3]2Ax[i4])
1293	dadg[i5][i5] = -Ax[i0]/denom
1294	for i0 in range(3):
1295	i1 = (i0+1)%3
1296	i2 = (i1+1)%3
1297	i3 = 5-i2
1298	for ij in range(6):
1299	dadg[i0][ij] = cell[i0](rcot[i2]dadg[i3][ij]/rsin[i2]-dsasdg[i1][ij]/ssin[i1])
1300	if ij == i0:
1301	dadg[i0][ij] = dadg[i0][ij]-0.5*cell[i0]/Ax[i0]
1302	dadg[i3][ij] = -dadg[i3][ij]rsin[2-i0]dpr
1303	sigMat = np.inner(dadg,np.inner(vcov,dadg.T))
1304	var = np.diag(sigMat)
1305	CS = np.where(var>0.,np.sqrt(var),0.)
1306	return [CS[0],CS[1],CS[2],CS[5],CS[4],CS[3],sigVol] #exchange sig(alp) & sig(gam) to get in right order
1307
1308	def SetPhaseData(parmDict,sigDict,Phases,RBIds,covData,RestraintDict=None,pFile=None):
1309	'needs a doc string'
1310
1311	def PrintAtomsAndSig(General,Atoms,atomsSig):
1312	print >>pFile,'\n Atoms:'
1313	line = ' name x y z frac Uiso U11 U22 U33 U12 U13 U23'
1314	if General['Type'] == 'magnetic':
1315	line += ' Mx My Mz'
1316	elif General['Type'] == 'macromolecular':
1317	line = ' res no residue chain '+line
1318	cx,ct,cs,cia = General['AtomPtrs']
1319	print >>pFile,line
1320	if General['Type'] == 'nuclear':
1321	print >>pFile,135*'-'
1322	fmt = {0:'%7s',1:'%7s',3:'%10.5f',4:'%10.5f',5:'%10.5f',6:'%8.3f',10:'%8.5f',
1323	11:'%8.5f',12:'%8.5f',13:'%8.5f',14:'%8.5f',15:'%8.5f',16:'%8.5f'}
1324	noFXsig = {3:[10' ','%10s'],4:[10' ','%10s'],5:[10' ','%10s'],6:[8' ','%8s']}
1325	for atyp in General['AtomTypes']: #zero composition
1326	General['NoAtoms'][atyp] = 0.
1327	for i,at in enumerate(Atoms):
1328	General['NoAtoms'][at[1]] += at[6]*float(at[8]) #new composition
1329	name = fmt[0]%(at[0])+fmt[1]%(at[1])+':'
1330	valstr = ' values:'
1331	sigstr = ' sig :'
1332	for ind in [3,4,5,6]:
1333	sigind = str(i)+':'+str(ind)
1334	valstr += fmt[ind]%(at[ind])
1335	if sigind in atomsSig:
1336	sigstr += fmt[ind]%(atomsSig[sigind])
1337	else:
1338	sigstr += noFXsig[ind][1]%(noFXsig[ind][0])
1339	if at[9] == 'I':
1340	valstr += fmt[10]%(at[10])
1341	if str(i)+':10' in atomsSig:
1342	sigstr += fmt[10]%(atomsSig[str(i)+':10'])
1343	else:
1344	sigstr += 8*' '
1345	else:
1346	valstr += 8*' '
1347	sigstr += 8*' '
1348	for ind in [11,12,13,14,15,16]:
1349	sigind = str(i)+':'+str(ind)
1350	valstr += fmt[ind]%(at[ind])
1351	if sigind in atomsSig:
1352	sigstr += fmt[ind]%(atomsSig[sigind])
1353	else:
1354	sigstr += 8*' '
1355	print >>pFile,name
1356	print >>pFile,valstr
1357	print >>pFile,sigstr
1358
1359	def PrintRBObjPOAndSig(rbfx,rbsx):
1360	namstr = ' names :'
1361	valstr = ' values:'
1362	sigstr = ' esds :'
1363	for i,px in enumerate(['Px:','Py:','Pz:']):
1364	name = pfx+rbfx+px+rbsx
1365	namstr += '%12s'%('Pos '+px[1])
1366	valstr += '%12.5f'%(parmDict[name])
1367	if name in sigDict:
1368	sigstr += '%12.5f'%(sigDict[name])
1369	else:
1370	sigstr += 12*' '
1371	for i,po in enumerate(['Oa:','Oi:','Oj:','Ok:']):
1372	name = pfx+rbfx+po+rbsx
1373	namstr += '%12s'%('Ori '+po[1])
1374	valstr += '%12.5f'%(parmDict[name])
1375	if name in sigDict:
1376	sigstr += '%12.5f'%(sigDict[name])
1377	else:
1378	sigstr += 12*' '
1379	print >>pFile,namstr
1380	print >>pFile,valstr
1381	print >>pFile,sigstr
1382
1383	def PrintRBObjTLSAndSig(rbfx,rbsx,TLS):
1384	namstr = ' names :'
1385	valstr = ' values:'
1386	sigstr = ' esds :'
1387	if 'N' not in TLS:
1388	print >>pFile,' Thermal motion:'
1389	if 'T' in TLS:
1390	for i,pt in enumerate(['T11:','T22:','T33:','T12:','T13:','T23:']):
1391	name = pfx+rbfx+pt+rbsx
1392	namstr += '%12s'%(pt[:3])
1393	valstr += '%12.5f'%(parmDict[name])
1394	if name in sigDict:
1395	sigstr += '%12.5f'%(sigDict[name])
1396	else:
1397	sigstr += 12*' '
1398	print >>pFile,namstr
1399	print >>pFile,valstr
1400	print >>pFile,sigstr
1401	if 'L' in TLS:
1402	namstr = ' names :'
1403	valstr = ' values:'
1404	sigstr = ' esds :'
1405	for i,pt in enumerate(['L11:','L22:','L33:','L12:','L13:','L23:']):
1406	name = pfx+rbfx+pt+rbsx
1407	namstr += '%12s'%(pt[:3])
1408	valstr += '%12.3f'%(parmDict[name])
1409	if name in sigDict:
1410	sigstr += '%12.3f'%(sigDict[name])
1411	else:
1412	sigstr += 12*' '
1413	print >>pFile,namstr
1414	print >>pFile,valstr
1415	print >>pFile,sigstr
1416	if 'S' in TLS:
1417	namstr = ' names :'
1418	valstr = ' values:'
1419	sigstr = ' esds :'
1420	for i,pt in enumerate(['S12:','S13:','S21:','S23:','S31:','S32:','SAA:','SBB:']):
1421	name = pfx+rbfx+pt+rbsx
1422	namstr += '%12s'%(pt[:3])
1423	valstr += '%12.4f'%(parmDict[name])
1424	if name in sigDict:
1425	sigstr += '%12.4f'%(sigDict[name])
1426	else:
1427	sigstr += 12*' '
1428	print >>pFile,namstr
1429	print >>pFile,valstr
1430	print >>pFile,sigstr
1431	if 'U' in TLS:
1432	name = pfx+rbfx+'U:'+rbsx
1433	namstr = ' names :'+'%12s'%('Uiso')
1434	valstr = ' values:'+'%12.5f'%(parmDict[name])
1435	if name in sigDict:
1436	sigstr = ' esds :'+'%12.5f'%(sigDict[name])
1437	else:
1438	sigstr = ' esds :'+12*' '
1439	print >>pFile,namstr
1440	print >>pFile,valstr
1441	print >>pFile,sigstr
1442
1443	def PrintRBObjTorAndSig(rbsx):
1444	namstr = ' names :'
1445	valstr = ' values:'
1446	sigstr = ' esds :'
1447	nTors = len(RBObj['Torsions'])
1448	if nTors:
1449	print >>pFile,' Torsions:'
1450	for it in range(nTors):
1451	name = pfx+'RBRTr;'+str(it)+':'+rbsx
1452	namstr += '%12s'%('Tor'+str(it))
1453	valstr += '%12.4f'%(parmDict[name])
1454	if name in sigDict:
1455	sigstr += '%12.4f'%(sigDict[name])
1456	print >>pFile,namstr
1457	print >>pFile,valstr
1458	print >>pFile,sigstr
1459
1460	def PrintSHtextureAndSig(textureData,SHtextureSig):
1461	print >>pFile,'\n Spherical harmonics texture: Order:' + str(textureData['Order'])
1462	names = ['omega','chi','phi']
1463	namstr = ' names :'
1464	ptstr = ' values:'
1465	sigstr = ' esds :'
1466	for name in names:
1467	namstr += '%12s'%(name)
1468	ptstr += '%12.3f'%(textureData['Sample '+name][1])
1469	if 'Sample '+name in SHtextureSig:
1470	sigstr += '%12.3f'%(SHtextureSig['Sample '+name])
1471	else:
1472	sigstr += 12*' '
1473	print >>pFile,namstr
1474	print >>pFile,ptstr
1475	print >>pFile,sigstr
1476	print >>pFile,'\n Texture coefficients:'
1477	namstr = ' names :'
1478	ptstr = ' values:'
1479	sigstr = ' esds :'
1480	SHcoeff = textureData['SH Coeff'][1]
1481	for name in SHcoeff:
1482	namstr += '%12s'%(name)
1483	ptstr += '%12.3f'%(SHcoeff[name])
1484	if name in SHtextureSig:
1485	sigstr += '%12.3f'%(SHtextureSig[name])
1486	else:
1487	sigstr += 12*' '
1488	print >>pFile,namstr
1489	print >>pFile,ptstr
1490	print >>pFile,sigstr
1491
1492	print >>pFile,'\n Phases:'
1493	for phase in Phases:
1494	print >>pFile,' Result for phase: ',phase
1495	Phase = Phases[phase]
1496	General = Phase['General']
1497	SGData = General['SGData']
1498	Atoms = Phase['Atoms']
1499	AtLookup = G2mth.FillAtomLookUp(Atoms)
1500	cell = General['Cell']
1501	pId = Phase['pId']
1502	pfx = str(pId)+'::'
1503	if cell[0]:
1504	A,sigA = cellFill(pfx,SGData,parmDict,sigDict)
1505	cellSig = getCellEsd(pfx,SGData,A,covData) #includes sigVol
1506	print >>pFile,' Reciprocal metric tensor: '
1507	ptfmt = "%15.9f"
1508	names = ['A11','A22','A33','A12','A13','A23']
1509	namstr = ' names :'
1510	ptstr = ' values:'
1511	sigstr = ' esds :'
1512	for name,a,siga in zip(names,A,sigA):
1513	namstr += '%15s'%(name)
1514	ptstr += ptfmt%(a)
1515	if siga:
1516	sigstr += ptfmt%(siga)
1517	else:
1518	sigstr += 15*' '
1519	print >>pFile,namstr
1520	print >>pFile,ptstr
1521	print >>pFile,sigstr
1522	cell[1:7] = G2lat.A2cell(A)
1523	cell[7] = G2lat.calc_V(A)
1524	print >>pFile,' New unit cell:'
1525	ptfmt = ["%12.6f","%12.6f","%12.6f","%12.4f","%12.4f","%12.4f","%12.3f"]
1526	names = ['a','b','c','alpha','beta','gamma','Volume']
1527	namstr = ' names :'
1528	ptstr = ' values:'
1529	sigstr = ' esds :'
1530	for name,fmt,a,siga in zip(names,ptfmt,cell[1:8],cellSig):
1531	namstr += '%12s'%(name)
1532	ptstr += fmt%(a)
1533	if siga:
1534	sigstr += fmt%(siga)
1535	else:
1536	sigstr += 12*' '
1537	print >>pFile,namstr
1538	print >>pFile,ptstr
1539	print >>pFile,sigstr
1540
1541	General['Mass'] = 0.
1542	if Phase['General'].get('doPawley'):
1543	pawleyRef = Phase['Pawley ref']
1544	for i,refl in enumerate(pawleyRef):
1545	key = pfx+'PWLref:'+str(i)
1546	refl[6] = parmDict[key]
1547	if key in sigDict:
1548	refl[7] = sigDict[key]
1549	else:
1550	refl[7] = 0
1551	else:
1552	VRBIds = RBIds['Vector']
1553	RRBIds = RBIds['Residue']
1554	RBModels = Phase['RBModels']
1555	for irb,RBObj in enumerate(RBModels.get('Vector',[])):
1556	jrb = VRBIds.index(RBObj['RBId'])
1557	rbsx = str(irb)+':'+str(jrb)
1558	print >>pFile,' Vector rigid body parameters:'
1559	PrintRBObjPOAndSig('RBV',rbsx)
1560	PrintRBObjTLSAndSig('RBV',rbsx,RBObj['ThermalMotion'][0])
1561	for irb,RBObj in enumerate(RBModels.get('Residue',[])):
1562	jrb = RRBIds.index(RBObj['RBId'])
1563	rbsx = str(irb)+':'+str(jrb)
1564	print >>pFile,' Residue rigid body parameters:'
1565	PrintRBObjPOAndSig('RBR',rbsx)
1566	PrintRBObjTLSAndSig('RBR',rbsx,RBObj['ThermalMotion'][0])
1567	PrintRBObjTorAndSig(rbsx)
1568	atomsSig = {}
1569	if General['Type'] == 'nuclear': #this needs macromolecular variant!
1570	for i,at in enumerate(Atoms):
1571	names = {3:pfx+'Ax:'+str(i),4:pfx+'Ay:'+str(i),5:pfx+'Az:'+str(i),6:pfx+'Afrac:'+str(i),
1572	10:pfx+'AUiso:'+str(i),11:pfx+'AU11:'+str(i),12:pfx+'AU22:'+str(i),13:pfx+'AU33:'+str(i),
1573	14:pfx+'AU12:'+str(i),15:pfx+'AU13:'+str(i),16:pfx+'AU23:'+str(i)}
1574	for ind in [3,4,5,6]:
1575	at[ind] = parmDict[names[ind]]
1576	if ind in [3,4,5]:
1577	name = names[ind].replace('A','dA')
1578	else:
1579	name = names[ind]
1580	if name in sigDict:
1581	atomsSig[str(i)+':'+str(ind)] = sigDict[name]
1582	if at[9] == 'I':
1583	at[10] = parmDict[names[10]]
1584	if names[10] in sigDict:
1585	atomsSig[str(i)+':10'] = sigDict[names[10]]
1586	else:
1587	for ind in [11,12,13,14,15,16]:
1588	at[ind] = parmDict[names[ind]]
1589	if names[ind] in sigDict:
1590	atomsSig[str(i)+':'+str(ind)] = sigDict[names[ind]]
1591	ind = General['AtomTypes'].index(at[1])
1592	General['Mass'] += General['AtomMass'][ind]at[6]at[8]
1593	PrintAtomsAndSig(General,Atoms,atomsSig)
1594
1595	textureData = General['SH Texture']
1596	if textureData['Order']:
1597	SHtextureSig = {}
1598	for name in ['omega','chi','phi']:
1599	aname = pfx+'SH '+name
1600	textureData['Sample '+name][1] = parmDict[aname]
1601	if aname in sigDict:
1602	SHtextureSig['Sample '+name] = sigDict[aname]
1603	for name in textureData['SH Coeff'][1]:
1604	aname = pfx+name
1605	textureData['SH Coeff'][1][name] = parmDict[aname]
1606	if aname in sigDict:
1607	SHtextureSig[name] = sigDict[aname]
1608	PrintSHtextureAndSig(textureData,SHtextureSig)
1609	if phase in RestraintDict:
1610	PrintRestraints(cell[1:7],SGData,General['AtomPtrs'],Atoms,AtLookup,
1611	textureData,RestraintDict[phase],pFile)
1612
1613	################################################################################
1614	##### Histogram & Phase data
1615	################################################################################
1616
1617	def GetHistogramPhaseData(Phases,Histograms,Print=True,pFile=None,resetRefList=True):
1618	'''Loads the HAP histogram/phase information into dicts
1619
1620	:param dict Phases: phase information
1621	:param dict Histograms: Histogram information
1622	:param bool Print: prints information as it is read
1623	:param file pFile: file object to print to (the default, None causes printing to the console)
1624	:param bool resetRefList: Should the contents of the Reflection List be initialized
1625	on loading. The default, True, initializes the Reflection List as it is loaded.
1626
1627	:returns: (hapVary,hapDict,controlDict)
1628	* hapVary: list of refined variables
1629	* hapDict: dict with refined variables and their values
1630	* controlDict: dict with computation controls (?)
1631	'''
1632
1633	def PrintSize(hapData):
1634	if hapData[0] in ['isotropic','uniaxial']:
1635	line = '\n Size model : %9s'%(hapData[0])
1636	line += ' equatorial:'+'%12.5f'%(hapData[1][0])+' Refine? '+str(hapData[2][0])
1637	if hapData[0] == 'uniaxial':
1638	line += ' axial:'+'%12.5f'%(hapData[1][1])+' Refine? '+str(hapData[2][1])
1639	line += '\n\t LG mixing coeff.: %12.4f'%(hapData[1][2])+' Refine? '+str(hapData[2][2])
1640	print >>pFile,line
1641	else:
1642	print >>pFile,'\n Size model : %s'%(hapData[0])+ \
1643	'\n\t LG mixing coeff.:%12.4f'%(hapData[1][2])+' Refine? '+str(hapData[2][2])
1644	Snames = ['S11','S22','S33','S12','S13','S23']
1645	ptlbls = ' names :'
1646	ptstr = ' values:'
1647	varstr = ' refine:'
1648	for i,name in enumerate(Snames):
1649	ptlbls += '%12s' % (name)
1650	ptstr += '%12.6f' % (hapData[4][i])
1651	varstr += '%12s' % (str(hapData[5][i]))
1652	print >>pFile,ptlbls
1653	print >>pFile,ptstr
1654	print >>pFile,varstr
1655
1656	def PrintMuStrain(hapData,SGData):
1657	if hapData[0] in ['isotropic','uniaxial']:
1658	line = '\n Mustrain model: %9s'%(hapData[0])
1659	line += ' equatorial:'+'%12.1f'%(hapData[1][0])+' Refine? '+str(hapData[2][0])
1660	if hapData[0] == 'uniaxial':
1661	line += ' axial:'+'%12.1f'%(hapData[1][1])+' Refine? '+str(hapData[2][1])
1662	line +='\n\t LG mixing coeff.:%12.4f'%(hapData[1][2])+' Refine? '+str(hapData[2][2])
1663	print >>pFile,line
1664	else:
1665	print >>pFile,'\n Mustrain model: %s'%(hapData[0])+ \
1666	'\n\t LG mixing coeff.:%12.4f'%(hapData[1][2])+' Refine? '+str(hapData[2][2])
1667	Snames = G2spc.MustrainNames(SGData)
1668	ptlbls = ' names :'
1669	ptstr = ' values:'
1670	varstr = ' refine:'
1671	for i,name in enumerate(Snames):
1672	ptlbls += '%12s' % (name)
1673	ptstr += '%12.6f' % (hapData[4][i])
1674	varstr += '%12s' % (str(hapData[5][i]))
1675	print >>pFile,ptlbls
1676	print >>pFile,ptstr
1677	print >>pFile,varstr
1678
1679	def PrintHStrain(hapData,SGData):
1680	print >>pFile,'\n Hydrostatic/elastic strain: '
1681	Hsnames = G2spc.HStrainNames(SGData)
1682	ptlbls = ' names :'
1683	ptstr = ' values:'
1684	varstr = ' refine:'
1685	for i,name in enumerate(Hsnames):
1686	ptlbls += '%12s' % (name)
1687	ptstr += '%12.4g' % (hapData[0][i])
1688	varstr += '%12s' % (str(hapData[1][i]))
1689	print >>pFile,ptlbls
1690	print >>pFile,ptstr
1691	print >>pFile,varstr
1692
1693	def PrintSHPO(hapData):
1694	print >>pFile,'\n Spherical harmonics preferred orientation: Order:' + \
1695	str(hapData[4])+' Refine? '+str(hapData[2])
1696	ptlbls = ' names :'
1697	ptstr = ' values:'
1698	for item in hapData[5]:
1699	ptlbls += '%12s'%(item)
1700	ptstr += '%12.3f'%(hapData[5][item])
1701	print >>pFile,ptlbls
1702	print >>pFile,ptstr
1703
1704	def PrintBabinet(hapData):
1705	print >>pFile,'\n Babinet form factor modification: '
1706	ptlbls = ' names :'
1707	ptstr = ' values:'
1708	varstr = ' refine:'
1709	for name in ['BabA','BabU']:
1710	ptlbls += '%12s' % (name)
1711	ptstr += '%12.6f' % (hapData[name][0])
1712	varstr += '%12s' % (str(hapData[name][1]))
1713	print >>pFile,ptlbls
1714	print >>pFile,ptstr
1715	print >>pFile,varstr
1716
1717
1718	hapDict = {}
1719	hapVary = []
1720	controlDict = {}
1721	poType = {}
1722	poAxes = {}
1723	spAxes = {}
1724	spType = {}
1725
1726	for phase in Phases:
1727	HistoPhase = Phases[phase]['Histograms']
1728	SGData = Phases[phase]['General']['SGData']
1729	cell = Phases[phase]['General']['Cell'][1:7]
1730	A = G2lat.cell2A(cell)
1731	pId = Phases[phase]['pId']
1732	histoList = HistoPhase.keys()
1733	histoList.sort()
1734	for histogram in histoList:
1735	try:
1736	Histogram = Histograms[histogram]
1737	except KeyError:
1738	#skip if histogram not included e.g. in a sequential refinement
1739	continue
1740	hapData = HistoPhase[histogram]
1741	hId = Histogram['hId']
1742	if 'PWDR' in histogram:
1743	limits = Histogram['Limits'][1]
1744	inst = Histogram['Instrument Parameters'][0]
1745	Zero = inst['Zero'][0]
1746	if 'C' in inst['Type'][1]:
1747	try:
1748	wave = inst['Lam'][1]
1749	except KeyError:
1750	wave = inst['Lam1'][1]
1751	dmin = wave/(2.0*sind(limits[1]/2.0))
1752	elif 'T' in inst['Type'][0]:
1753	dmin = limits[0]/inst['difC'][1]
1754	pfx = str(pId)+':'+str(hId)+':'
1755	for item in ['Scale','Extinction']:
1756	hapDict[pfx+item] = hapData[item][0]
1757	if hapData[item][1]:
1758	hapVary.append(pfx+item)
1759	names = G2spc.HStrainNames(SGData)
1760	HSvals = []
1761	for i,name in enumerate(names):
1762	hapDict[pfx+name] = hapData['HStrain'][0][i]
1763	HSvals.append(hapDict[pfx+name])
1764	if hapData['HStrain'][1][i]:
1765	hapVary.append(pfx+name)
1766	DIJS = G2spc.HStrainVals(HSvals,SGData)
1767	controlDict[pfx+'poType'] = hapData['Pref.Ori.'][0]
1768	if hapData['Pref.Ori.'][0] == 'MD':
1769	hapDict[pfx+'MD'] = hapData['Pref.Ori.'][1]
1770	controlDict[pfx+'MDAxis'] = hapData['Pref.Ori.'][3]
1771	if hapData['Pref.Ori.'][2]:
1772	hapVary.append(pfx+'MD')
1773	else: #'SH' spherical harmonics
1774	controlDict[pfx+'SHord'] = hapData['Pref.Ori.'][4]
1775	controlDict[pfx+'SHncof'] = len(hapData['Pref.Ori.'][5])
1776	for item in hapData['Pref.Ori.'][5]:
1777	hapDict[pfx+item] = hapData['Pref.Ori.'][5][item]
1778	if hapData['Pref.Ori.'][2]:
1779	hapVary.append(pfx+item)
1780	for item in ['Mustrain','Size']:
1781	controlDict[pfx+item+'Type'] = hapData[item][0]
1782	hapDict[pfx+item+';mx'] = hapData[item][1][2]
1783	if hapData[item][2][2]:
1784	hapVary.append(pfx+item+';mx')
1785	if hapData[item][0] in ['isotropic','uniaxial']:
1786	hapDict[pfx+item+';i'] = hapData[item][1][0]
1787	if hapData[item][2][0]:
1788	hapVary.append(pfx+item+';i')
1789	if hapData[item][0] == 'uniaxial':
1790	controlDict[pfx+item+'Axis'] = hapData[item][3]
1791	hapDict[pfx+item+';a'] = hapData[item][1][1]
1792	if hapData[item][2][1]:
1793	hapVary.append(pfx+item+';a')
1794	else: #generalized for mustrain or ellipsoidal for size
1795	Nterms = len(hapData[item][4])
1796	if item == 'Mustrain':
1797	names = G2spc.MustrainNames(SGData)
1798	pwrs = []
1799	for name in names:
1800	h,k,l = name[1:]
1801	pwrs.append([int(h),int(k),int(l)])
1802	controlDict[pfx+'MuPwrs'] = pwrs
1803	for i in range(Nterms):
1804	sfx = ':'+str(i)
1805	hapDict[pfx+item+sfx] = hapData[item][4][i]
1806	if hapData[item][5][i]:
1807	hapVary.append(pfx+item+sfx)
1808	for bab in ['BabA','BabU']:
1809	hapDict[pfx+bab] = hapData['Babinet'][bab][0]
1810	if hapData['Babinet'][bab][1]:
1811	hapVary.append(pfx+bab)
1812
1813	if Print:
1814	print >>pFile,'\n Phase: ',phase,' in histogram: ',histogram
1815	print >>pFile,135*'-'
1816	print >>pFile,' Phase fraction : %10.4f'%(hapData['Scale'][0]),' Refine?',hapData['Scale'][1]
1817	print >>pFile,' Extinction coeff: %10.4f'%(hapData['Extinction'][0]),' Refine?',hapData['Extinction'][1]
1818	if hapData['Pref.Ori.'][0] == 'MD':
1819	Ax = hapData['Pref.Ori.'][3]
1820	print >>pFile,' March-Dollase PO: %10.4f'%(hapData['Pref.Ori.'][1]),' Refine?',hapData['Pref.Ori.'][2], \
1821	' Axis: %d %d %d'%(Ax[0],Ax[1],Ax[2])
1822	else: #'SH' for spherical harmonics
1823	PrintSHPO(hapData['Pref.Ori.'])
1824	PrintSize(hapData['Size'])
1825	PrintMuStrain(hapData['Mustrain'],SGData)
1826	PrintHStrain(hapData['HStrain'],SGData)
1827	if hapData['Babinet']['BabA'][0]:
1828	PrintBabinet(hapData['Babinet'])
1829	HKLd = np.array(G2lat.GenHLaue(dmin,SGData,A)) #+DIJS
1830	if resetRefList:
1831	refList = []
1832	Uniq = []
1833	Phi = []
1834	for h,k,l,d in HKLd:
1835	ext,mul,uniq,phi = G2spc.GenHKLf([h,k,l],SGData)
1836	mul *= 2 # for powder overlap of Friedel pairs
1837	if ext:
1838	continue
1839	if 'C' in inst['Type'][0]:
1840	pos = 2.0asind(wave/(2.0d))+Zero
1841	if limits[0] < pos < limits[1]:
1842	refList.append([h,k,l,mul,d, pos,0.0,0.0,0.0,0.0, 0.0,0.0,1.0,1.0,1.0])
1843	#... sig,gam,fotsq,fctsq, phase,icorr,prfo,abs,ext
1844	Uniq.append(uniq)
1845	Phi.append(phi)
1846	elif 'T' in inst['Type'][0]:
1847	pos = inst['difC'][1]d+inst['difA'][1]d**2+inst['difB'][1]/d+Zero
1848	if limits[0] < pos < limits[1]:
1849	wave = inst['difC'][1]d/(252.816inst['fltPath'][0])
1850	refList.append([h,k,l,mul,d, pos,0.0,0.0,0.0,0.0, 0.0,0.0,0.0,0.0,wave, 1.0,1.0,1.0])
1851	# ... sig,gam,fotsq,fctsq, phase,icorr,alp,bet,wave, prfo,abs,ext
1852	Uniq.append(uniq)
1853	Phi.append(phi)
1854	Histogram['Reflection Lists'][phase] = {'RefList':np.array(refList),'FF':{},'Type':inst['Type'][0]}
1855	elif 'HKLF' in histogram:
1856	inst = Histogram['Instrument Parameters'][0]
1857	hId = Histogram['hId']
1858	hfx = ':%d:'%(hId)
1859	for item in inst:
1860	if type(inst) is not list and item != 'Type': continue # skip over non-refined items (such as InstName)
1861	hapDict[hfx+item] = inst[item][1]
1862	pfx = str(pId)+':'+str(hId)+':'
1863	hapDict[pfx+'Scale'] = hapData['Scale'][0]
1864	if hapData['Scale'][1]:
1865	hapVary.append(pfx+'Scale')
1866
1867	extApprox,extType,extParms = hapData['Extinction']
1868	controlDict[pfx+'EType'] = extType
1869	controlDict[pfx+'EApprox'] = extApprox
1870	if 'C' in inst['Type'][0]:
1871	controlDict[pfx+'Tbar'] = extParms['Tbar']
1872	controlDict[pfx+'Cos2TM'] = extParms['Cos2TM']
1873	if 'Primary' in extType:
1874	Ekey = ['Ep',]
1875	elif 'I & II' in extType:
1876	Ekey = ['Eg','Es']
1877	elif 'Secondary Type II' == extType:
1878	Ekey = ['Es',]
1879	elif 'Secondary Type I' == extType:
1880	Ekey = ['Eg',]
1881	else: #'None'
1882	Ekey = []
1883	for eKey in Ekey:
1884	hapDict[pfx+eKey] = extParms[eKey][0]
1885	if extParms[eKey][1]:
1886	hapVary.append(pfx+eKey)
1887	for bab in ['BabA','BabU']:
1888	hapDict[pfx+bab] = hapData['Babinet'][bab][0]
1889	if hapData['Babinet'][bab][1]:
1890	hapVary.append(pfx+bab)
1891	if Print:
1892	print >>pFile,'\n Phase: ',phase,' in histogram: ',histogram
1893	print >>pFile,135*'-'
1894	print >>pFile,' Scale factor : %10.4f'%(hapData['Scale'][0]),' Refine?',hapData['Scale'][1]
1895	if extType != 'None':
1896	print >>pFile,' Extinction Type: %15s'%(extType),' approx: %10s'%(extApprox)
1897	text = ' Parameters :'
1898	for eKey in Ekey:
1899	text += ' %4s : %10.3e Refine? '%(eKey,extParms[eKey][0])+str(extParms[eKey][1])
1900	print >>pFile,text
1901	if hapData['Babinet']['BabA'][0]:
1902	PrintBabinet(hapData['Babinet'])
1903	Histogram['Reflection Lists'] = phase
1904
1905	return hapVary,hapDict,controlDict
1906
1907	def SetHistogramPhaseData(parmDict,sigDict,Phases,Histograms,Print=True,pFile=None):
1908	'needs a doc string'
1909
1910	def PrintSizeAndSig(hapData,sizeSig):
1911	line = '\n Size model: %9s'%(hapData[0])
1912	refine = False
1913	if hapData[0] in ['isotropic','uniaxial']:
1914	line += ' equatorial:%12.5f'%(hapData[1][0])
1915	if sizeSig[0][0]:
1916	line += ', sig:%8.4f'%(sizeSig[0][0])
1917	refine = True
1918	if hapData[0] == 'uniaxial':
1919	line += ' axial:%12.4f'%(hapData[1][1])
1920	if sizeSig[0][1]:
1921	refine = True
1922	line += ', sig:%8.4f'%(sizeSig[0][1])
1923	line += ' LG mix coeff.:%12.4f'%(hapData[1][2])
1924	if sizeSig[0][2]:
1925	refine = True
1926	line += ', sig:%8.4f'%(sizeSig[0][2])
1927	if refine:
1928	print >>pFile,line
1929	else:
1930	line += ' LG mix coeff.:%12.4f'%(hapData[1][2])
1931	if sizeSig[0][2]:
1932	refine = True
1933	line += ', sig:%8.4f'%(sizeSig[0][2])
1934	Snames = ['S11','S22','S33','S12','S13','S23']
1935	ptlbls = ' name :'
1936	ptstr = ' value :'
1937	sigstr = ' sig :'
1938	for i,name in enumerate(Snames):
1939	ptlbls += '%12s' % (name)
1940	ptstr += '%12.6f' % (hapData[4][i])
1941	if sizeSig[1][i]:
1942	refine = True
1943	sigstr += '%12.6f' % (sizeSig[1][i])
1944	else:
1945	sigstr += 12*' '
1946	if refine:
1947	print >>pFile,line
1948	print >>pFile,ptlbls
1949	print >>pFile,ptstr
1950	print >>pFile,sigstr
1951
1952	def PrintMuStrainAndSig(hapData,mustrainSig,SGData):
1953	line = '\n Mustrain model: %9s'%(hapData[0])
1954	refine = False
1955	if hapData[0] in ['isotropic','uniaxial']:
1956	line += ' equatorial:%12.1f'%(hapData[1][0])
1957	if mustrainSig[0][0]:
1958	line += ', sig:%8.1f'%(mustrainSig[0][0])
1959	refine = True
1960	if hapData[0] == 'uniaxial':
1961	line += ' axial:%12.1f'%(hapData[1][1])
1962	if mustrainSig[0][1]:
1963	line += ', sig:%8.1f'%(mustrainSig[0][1])
1964	line += ' LG mix coeff.:%12.4f'%(hapData[1][2])
1965	if mustrainSig[0][2]:
1966	refine = True
1967	line += ', sig:%8.3f'%(mustrainSig[0][2])
1968	if refine:
1969	print >>pFile,line
1970	else:
1971	line += ' LG mix coeff.:%12.4f'%(hapData[1][2])
1972	if mustrainSig[0][2]:
1973	refine = True
1974	line += ', sig:%8.3f'%(mustrainSig[0][2])
1975	Snames = G2spc.MustrainNames(SGData)
1976	ptlbls = ' name :'
1977	ptstr = ' value :'
1978	sigstr = ' sig :'
1979	for i,name in enumerate(Snames):
1980	ptlbls += '%12s' % (name)
1981	ptstr += '%12.6f' % (hapData[4][i])
1982	if mustrainSig[1][i]:
1983	refine = True
1984	sigstr += '%12.6f' % (mustrainSig[1][i])
1985	else:
1986	sigstr += 12*' '
1987	if refine:
1988	print >>pFile,line
1989	print >>pFile,ptlbls
1990	print >>pFile,ptstr
1991	print >>pFile,sigstr
1992
1993	def PrintHStrainAndSig(hapData,strainSig,SGData):
1994	Hsnames = G2spc.HStrainNames(SGData)
1995	ptlbls = ' name :'
1996	ptstr = ' value :'
1997	sigstr = ' sig :'
1998	refine = False
1999	for i,name in enumerate(Hsnames):
2000	ptlbls += '%12s' % (name)
2001	ptstr += '%12.4g' % (hapData[0][i])
2002	if name in strainSig:
2003	refine = True
2004	sigstr += '%12.4g' % (strainSig[name])
2005	else:
2006	sigstr += 12*' '
2007	if refine:
2008	print >>pFile,'\n Hydrostatic/elastic strain: '
2009	print >>pFile,ptlbls
2010	print >>pFile,ptstr
2011	print >>pFile,sigstr
2012
2013	def PrintSHPOAndSig(pfx,hapData,POsig):
2014	print >>pFile,'\n Spherical harmonics preferred orientation: Order:'+str(hapData[4])
2015	ptlbls = ' names :'
2016	ptstr = ' values:'
2017	sigstr = ' sig :'
2018	for item in hapData[5]:
2019	ptlbls += '%12s'%(item)
2020	ptstr += '%12.3f'%(hapData[5][item])
2021	if pfx+item in POsig:
2022	sigstr += '%12.3f'%(POsig[pfx+item])
2023	else:
2024	sigstr += 12*' '
2025	print >>pFile,ptlbls
2026	print >>pFile,ptstr
2027	print >>pFile,sigstr
2028
2029	def PrintExtAndSig(pfx,hapData,ScalExtSig):
2030	print >>pFile,'\n Single crystal extinction: Type: ',hapData[0],' Approx: ',hapData[1]
2031	text = ''
2032	for item in hapData[2]:
2033	if pfx+item in ScalExtSig:
2034	text += ' %s: '%(item)
2035	text += '%12.2e'%(hapData[2][item][0])
2036	if pfx+item in ScalExtSig:
2037	text += ' sig: %12.2e'%(ScalExtSig[pfx+item])
2038	print >>pFile,text
2039
2040	def PrintBabinetAndSig(pfx,hapData,BabSig):
2041	print >>pFile,'\n Babinet form factor modification: '
2042	ptlbls = ' names :'
2043	ptstr = ' values:'
2044	sigstr = ' sig :'
2045	for item in hapData:
2046	ptlbls += '%12s'%(item)
2047	ptstr += '%12.3f'%(hapData[item][0])
2048	if pfx+item in BabSig:
2049	sigstr += '%12.3f'%(BabSig[pfx+item])
2050	else:
2051	sigstr += 12*' '
2052	print >>pFile,ptlbls
2053	print >>pFile,ptstr
2054	print >>pFile,sigstr
2055
2056	PhFrExtPOSig = {}
2057	SizeMuStrSig = {}
2058	ScalExtSig = {}
2059	BabSig = {}
2060	wtFrSum = {}
2061	for phase in Phases:
2062	HistoPhase = Phases[phase]['Histograms']
2063	General = Phases[phase]['General']
2064	SGData = General['SGData']
2065	pId = Phases[phase]['pId']
2066	histoList = HistoPhase.keys()
2067	histoList.sort()
2068	for histogram in histoList:
2069	try:
2070	Histogram = Histograms[histogram]
2071	except KeyError:
2072	#skip if histogram not included e.g. in a sequential refinement
2073	continue
2074	hapData = HistoPhase[histogram]
2075	hId = Histogram['hId']
2076	pfx = str(pId)+':'+str(hId)+':'
2077	if hId not in wtFrSum:
2078	wtFrSum[hId] = 0.
2079	if 'PWDR' in histogram:
2080	for item in ['Scale','Extinction']:
2081	hapData[item][0] = parmDict[pfx+item]
2082	if pfx+item in sigDict:
2083	PhFrExtPOSig.update({pfx+item:sigDict[pfx+item],})
2084	wtFrSum[hId] += hapData['Scale'][0]*General['Mass']
2085	if hapData['Pref.Ori.'][0] == 'MD':
2086	hapData['Pref.Ori.'][1] = parmDict[pfx+'MD']
2087	if pfx+'MD' in sigDict:
2088	PhFrExtPOSig.update({pfx+'MD':sigDict[pfx+'MD'],})
2089	else: #'SH' spherical harmonics
2090	for item in hapData['Pref.Ori.'][5]:
2091	hapData['Pref.Ori.'][5][item] = parmDict[pfx+item]
2092	if pfx+item in sigDict:
2093	PhFrExtPOSig.update({pfx+item:sigDict[pfx+item],})
2094	SizeMuStrSig.update({pfx+'Mustrain':[[0,0,0],[0 for i in range(len(hapData['Mustrain'][4]))]],
2095	pfx+'Size':[[0,0,0],[0 for i in range(len(hapData['Size'][4]))]],
2096	pfx+'HStrain':{}})
2097	for item in ['Mustrain','Size']:
2098	hapData[item][1][2] = parmDict[pfx+item+';mx']
2099	hapData[item][1][2] = min(1.,max(0.,hapData[item][1][2]))
2100	if pfx+item+';mx' in sigDict:
2101	SizeMuStrSig[pfx+item][0][2] = sigDict[pfx+item+';mx']
2102	if hapData[item][0] in ['isotropic','uniaxial']:
2103	hapData[item][1][0] = parmDict[pfx+item+';i']
2104	if item == 'Size':
2105	hapData[item][1][0] = min(10.,max(0.001,hapData[item][1][0]))
2106	if pfx+item+';i' in sigDict:
2107	SizeMuStrSig[pfx+item][0][0] = sigDict[pfx+item+';i']
2108	if hapData[item][0] == 'uniaxial':
2109	hapData[item][1][1] = parmDict[pfx+item+';a']
2110	if item == 'Size':
2111	hapData[item][1][1] = min(10.,max(0.001,hapData[item][1][1]))
2112	if pfx+item+';a' in sigDict:
2113	SizeMuStrSig[pfx+item][0][1] = sigDict[pfx+item+';a']
2114	else: #generalized for mustrain or ellipsoidal for size
2115	Nterms = len(hapData[item][4])
2116	for i in range(Nterms):
2117	sfx = ':'+str(i)
2118	hapData[item][4][i] = parmDict[pfx+item+sfx]
2119	if pfx+item+sfx in sigDict:
2120	SizeMuStrSig[pfx+item][1][i] = sigDict[pfx+item+sfx]
2121	names = G2spc.HStrainNames(SGData)
2122	for i,name in enumerate(names):
2123	hapData['HStrain'][0][i] = parmDict[pfx+name]
2124	if pfx+name in sigDict:
2125	SizeMuStrSig[pfx+'HStrain'][name] = sigDict[pfx+name]
2126	for name in ['BabA','BabU']:
2127	hapData['Babinet'][name][0] = parmDict[pfx+name]
2128	if pfx+name in sigDict:
2129	BabSig[pfx+name] = sigDict[pfx+name]
2130
2131	elif 'HKLF' in histogram:
2132	for item in ['Scale',]:
2133	if parmDict.get(pfx+item):
2134	hapData[item][0] = parmDict[pfx+item]
2135	if pfx+item in sigDict:
2136	ScalExtSig[pfx+item] = sigDict[pfx+item]
2137	for item in ['Ep','Eg','Es']:
2138	if parmDict.get(pfx+item):
2139	hapData['Extinction'][2][item][0] = parmDict[pfx+item]
2140	if pfx+item in sigDict:
2141	ScalExtSig[pfx+item] = sigDict[pfx+item]
2142	for name in ['BabA','BabU']:
2143	hapData['Babinet'][name][0] = parmDict[pfx+name]
2144	if pfx+name in sigDict:
2145	BabSig[pfx+name] = sigDict[pfx+name]
2146
2147	if Print:
2148	for phase in Phases:
2149	HistoPhase = Phases[phase]['Histograms']
2150	General = Phases[phase]['General']
2151	SGData = General['SGData']
2152	pId = Phases[phase]['pId']
2153	histoList = HistoPhase.keys()
2154	histoList.sort()
2155	for histogram in histoList:
2156	try:
2157	Histogram = Histograms[histogram]
2158	except KeyError:
2159	#skip if histogram not included e.g. in a sequential refinement
2160	continue
2161	print >>pFile,'\n Phase: ',phase,' in histogram: ',histogram
2162	print >>pFile,130*'-'
2163	hapData = HistoPhase[histogram]
2164	hId = Histogram['hId']
2165	Histogram['Residuals'][str(pId)+'::Name'] = phase
2166	pfx = str(pId)+':'+str(hId)+':'
2167	if 'PWDR' in histogram:
2168	print >>pFile,' Final refinement RF, RF^2 = %.2f%%, %.2f%% on %d reflections' \
2169	%(Histogram['Residuals'][pfx+'Rf'],Histogram['Residuals'][pfx+'Rf^2'],Histogram['Residuals'][pfx+'Nref'])
2170
2171	if pfx+'Scale' in PhFrExtPOSig:
2172	wtFr = hapData['Scale'][0]*General['Mass']/wtFrSum[hId]
2173	sigwtFr = PhFrExtPOSig[pfx+'Scale']*wtFr/hapData['Scale'][0]
2174	print >>pFile,' Phase fraction : %10.5f, sig %10.5f Weight fraction : %8.5f, sig %10.5f' \
2175	%(hapData['Scale'][0],PhFrExtPOSig[pfx+'Scale'],wtFr,sigwtFr)
2176	if pfx+'Extinction' in PhFrExtPOSig:
2177	print >>pFile,' Extinction coeff: %10.4f, sig %10.4f'%(hapData['Extinction'][0],PhFrExtPOSig[pfx+'Extinction'])
2178	if hapData['Pref.Ori.'][0] == 'MD':
2179	if pfx+'MD' in PhFrExtPOSig:
2180	print >>pFile,' March-Dollase PO: %10.4f, sig %10.4f'%(hapData['Pref.Ori.'][1],PhFrExtPOSig[pfx+'MD'])
2181	else:
2182	PrintSHPOAndSig(pfx,hapData['Pref.Ori.'],PhFrExtPOSig)
2183	PrintSizeAndSig(hapData['Size'],SizeMuStrSig[pfx+'Size'])
2184	PrintMuStrainAndSig(hapData['Mustrain'],SizeMuStrSig[pfx+'Mustrain'],SGData)
2185	PrintHStrainAndSig(hapData['HStrain'],SizeMuStrSig[pfx+'HStrain'],SGData)
2186	if len(BabSig):
2187	PrintBabinetAndSig(pfx,hapData['Babinet'],BabSig)
2188
2189	elif 'HKLF' in histogram:
2190	print >>pFile,' Final refinement RF, RF^2 = %.2f%%, %.2f%% on %d reflections' \
2191	%(Histogram['Residuals'][pfx+'Rf'],Histogram['Residuals'][pfx+'Rf^2'],Histogram['Residuals'][pfx+'Nref'])
2192	print >>pFile,' HKLF histogram weight factor = ','%.3f'%(Histogram['wtFactor'])
2193	if pfx+'Scale' in ScalExtSig:
2194	print >>pFile,' Scale factor : %10.4f, sig %10.4f'%(hapData['Scale'][0],ScalExtSig[pfx+'Scale'])
2195	if hapData['Extinction'][0] != 'None':
2196	PrintExtAndSig(pfx,hapData['Extinction'],ScalExtSig)
2197	if len(BabSig):
2198	PrintBabinetAndSig(pfx,hapData['Babinet'],BabSig)
2199
2200	################################################################################
2201	##### Histogram data
2202	################################################################################
2203
2204	def GetHistogramData(Histograms,Print=True,pFile=None):
2205	'needs a doc string'
2206
2207	def GetBackgroundParms(hId,Background):
2208	Back = Background[0]
2209	DebyePeaks = Background[1]
2210	bakType,bakFlag = Back[:2]
2211	backVals = Back[3:]
2212	backNames = [':'+str(hId)+':Back;'+str(i) for i in range(len(backVals))]
2213	backDict = dict(zip(backNames,backVals))
2214	backVary = []
2215	if bakFlag:
2216	backVary = backNames
2217	backDict[':'+str(hId)+':nDebye'] = DebyePeaks['nDebye']
2218	backDict[':'+str(hId)+':nPeaks'] = DebyePeaks['nPeaks']
2219	debyeDict = {}
2220	debyeList = []
2221	for i in range(DebyePeaks['nDebye']):
2222	debyeNames = [':'+str(hId)+':DebyeA;'+str(i),':'+str(hId)+':DebyeR;'+str(i),':'+str(hId)+':DebyeU;'+str(i)]
2223	debyeDict.update(dict(zip(debyeNames,DebyePeaks['debyeTerms'][i][::2])))
2224	debyeList += zip(debyeNames,DebyePeaks['debyeTerms'][i][1::2])
2225	debyeVary = []
2226	for item in debyeList:
2227	if item[1]:
2228	debyeVary.append(item[0])
2229	backDict.update(debyeDict)
2230	backVary += debyeVary
2231	peakDict = {}
2232	peakList = []
2233	for i in range(DebyePeaks['nPeaks']):
2234	peakNames = [':'+str(hId)+':BkPkpos;'+str(i),':'+str(hId)+ \
2235	':BkPkint;'+str(i),':'+str(hId)+':BkPksig;'+str(i),':'+str(hId)+':BkPkgam;'+str(i)]
2236	peakDict.update(dict(zip(peakNames,DebyePeaks['peaksList'][i][::2])))
2237	peakList += zip(peakNames,DebyePeaks['peaksList'][i][1::2])
2238	peakVary = []
2239	for item in peakList:
2240	if item[1]:
2241	peakVary.append(item[0])
2242	backDict.update(peakDict)
2243	backVary += peakVary
2244	return bakType,backDict,backVary
2245
2246	def GetInstParms(hId,Inst):
2247	dataType = Inst['Type'][0]
2248	instDict = {}
2249	insVary = []
2250	pfx = ':'+str(hId)+':'
2251	insKeys = Inst.keys()
2252	insKeys.sort()
2253	for item in insKeys:
2254	insName = pfx+item
2255	instDict[insName] = Inst[item][1]
2256	if len(Inst[item]) > 2 and Inst[item][2]:
2257	insVary.append(insName)
2258	if 'C' in dataType:
2259	instDict[pfx+'SH/L'] = max(instDict[pfx+'SH/L'],0.0005)
2260	return dataType,instDict,insVary
2261
2262	def GetSampleParms(hId,Sample):
2263	sampVary = []
2264	hfx = ':'+str(hId)+':'
2265	sampDict = {hfx+'Gonio. radius':Sample['Gonio. radius'],hfx+'Omega':Sample['Omega'],
2266	hfx+'Chi':Sample['Chi'],hfx+'Phi':Sample['Phi']}
2267	for key in ('Temperature','Pressure','FreePrm1','FreePrm2','FreePrm3'):
2268	if key in Sample:
2269	sampDict[hfx+key] = Sample[key]
2270	Type = Sample['Type']
2271	if 'Bragg' in Type: #Bragg-Brentano
2272	for item in ['Scale','Shift','Transparency','SurfRoughA','SurfRoughB']:
2273	sampDict[hfx+item] = Sample[item][0]
2274	if Sample[item][1]:
2275	sampVary.append(hfx+item)
2276	elif 'Debye' in Type: #Debye-Scherrer
2277	for item in ['Scale','Absorption','DisplaceX','DisplaceY']:
2278	sampDict[hfx+item] = Sample[item][0]
2279	if Sample[item][1]:
2280	sampVary.append(hfx+item)
2281	return Type,sampDict,sampVary
2282
2283	def PrintBackground(Background):
2284	Back = Background[0]
2285	DebyePeaks = Background[1]
2286	print >>pFile,'\n Background function: ',Back[0],' Refine?',bool(Back[1])
2287	line = ' Coefficients: '
2288	for i,back in enumerate(Back[3:]):
2289	line += '%10.3f'%(back)
2290	if i and not i%10:
2291	line += '\n'+15*' '
2292	print >>pFile,line
2293	if DebyePeaks['nDebye']:
2294	print >>pFile,'\n Debye diffuse scattering coefficients'
2295	parms = ['DebyeA','DebyeR','DebyeU']
2296	line = ' names : '
2297	for parm in parms:
2298	line += '%8s refine?'%(parm)
2299	print >>pFile,line
2300	for j,term in enumerate(DebyePeaks['debyeTerms']):
2301	line = ' term'+'%2d'%(j)+':'
2302	for i in range(3):
2303	line += '%10.3f %5s'%(term[2i],bool(term[2i+1]))
2304	print >>pFile,line
2305	if DebyePeaks['nPeaks']:
2306	print >>pFile,'\n Single peak coefficients'
2307	parms = ['BkPkpos','BkPkint','BkPksig','BkPkgam']
2308	line = ' names : '
2309	for parm in parms:
2310	line += '%8s refine?'%(parm)
2311	print >>pFile,line
2312	for j,term in enumerate(DebyePeaks['peaksList']):
2313	line = ' peak'+'%2d'%(j)+':'
2314	for i in range(4):
2315	line += '%10.3f %5s'%(term[2i],bool(term[2i+1]))
2316	print >>pFile,line
2317
2318	def PrintInstParms(Inst):
2319	print >>pFile,'\n Instrument Parameters:'
2320	insKeys = Inst.keys()
2321	insKeys.sort()
2322	iBeg = 0
2323	Ok = True
2324	while Ok:
2325	ptlbls = ' name :'
2326	ptstr = ' value :'
2327	varstr = ' refine:'
2328	iFin = min(iBeg+9,len(insKeys))
2329	for item in insKeys[iBeg:iFin]:
2330	if item not in ['Type','Source']:
2331	ptlbls += '%12s' % (item)
2332	ptstr += '%12.6f' % (Inst[item][1])
2333	if item in ['Lam1','Lam2','Azimuth','fltPath','2-theta',]:
2334	varstr += 12*' '
2335	else:
2336	varstr += '%12s' % (str(bool(Inst[item][2])))
2337	print >>pFile,ptlbls
2338	print >>pFile,ptstr
2339	print >>pFile,varstr
2340	iBeg = iFin
2341	if iBeg == len(insKeys):
2342	Ok = False
2343	else:
2344	print >>pFile,'\n'
2345
2346	def PrintSampleParms(Sample):
2347	print >>pFile,'\n Sample Parameters:'
2348	print >>pFile,' Goniometer omega = %.2f, chi = %.2f, phi = %.2f'% \
2349	(Sample['Omega'],Sample['Chi'],Sample['Phi'])
2350	ptlbls = ' name :'
2351	ptstr = ' value :'
2352	varstr = ' refine:'
2353	if 'Bragg' in Sample['Type']:
2354	for item in ['Scale','Shift','Transparency','SurfRoughA','SurfRoughB']:
2355	ptlbls += '%14s'%(item)
2356	ptstr += '%14.4f'%(Sample[item][0])
2357	varstr += '%14s'%(str(bool(Sample[item][1])))
2358
2359	elif 'Debye' in Type: #Debye-Scherrer
2360	for item in ['Scale','Absorption','DisplaceX','DisplaceY']:
2361	ptlbls += '%14s'%(item)
2362	ptstr += '%14.4f'%(Sample[item][0])
2363	varstr += '%14s'%(str(bool(Sample[item][1])))
2364
2365	print >>pFile,ptlbls
2366	print >>pFile,ptstr
2367	print >>pFile,varstr
2368
2369	histDict = {}
2370	histVary = []
2371	controlDict = {}
2372	histoList = Histograms.keys()
2373	histoList.sort()
2374	for histogram in histoList:
2375	if 'PWDR' in histogram:
2376	Histogram = Histograms[histogram]
2377	hId = Histogram['hId']
2378	pfx = ':'+str(hId)+':'
2379	controlDict[pfx+'wtFactor'] = Histogram['wtFactor']
2380	controlDict[pfx+'Limits'] = Histogram['Limits'][1]
2381	controlDict[pfx+'Exclude'] = Histogram['Limits'][2:]
2382	for excl in controlDict[pfx+'Exclude']:
2383	Histogram['Data'][0] = ma.masked_inside(Histogram['Data'][0],excl[0],excl[1])
2384	if controlDict[pfx+'Exclude']:
2385	ma.mask_rows(Histogram['Data'])
2386	Background = Histogram['Background']
2387	Type,bakDict,bakVary = GetBackgroundParms(hId,Background)
2388	controlDict[pfx+'bakType'] = Type
2389	histDict.update(bakDict)
2390	histVary += bakVary
2391
2392	Inst = Histogram['Instrument Parameters'][0]
2393	Type,instDict,insVary = GetInstParms(hId,Inst)
2394	controlDict[pfx+'histType'] = Type
2395	if 'XC' in Type:
2396	if pfx+'Lam1' in instDict:
2397	controlDict[pfx+'keV'] = 12.397639/instDict[pfx+'Lam1']
2398	else:
2399	controlDict[pfx+'keV'] = 12.397639/instDict[pfx+'Lam']
2400	histDict.update(instDict)
2401	histVary += insVary
2402
2403	Sample = Histogram['Sample Parameters']
2404	Type,sampDict,sampVary = GetSampleParms(hId,Sample)
2405	controlDict[pfx+'instType'] = Type
2406	histDict.update(sampDict)
2407	histVary += sampVary
2408
2409
2410	if Print:
2411	print >>pFile,'\n Histogram: ',histogram,' histogram Id: ',hId
2412	print >>pFile,135*'-'
2413	Units = {'C':' deg','T':' msec'}
2414	units = Units[controlDict[pfx+'histType'][2]]
2415	Limits = controlDict[pfx+'Limits']
2416	print >>pFile,' Instrument type: ',Sample['Type']
2417	print >>pFile,' Histogram limits: %8.2f%s to %8.2f%s'%(Limits[0],units,Limits[1],units)
2418	if len(controlDict[pfx+'Exclude']):
2419	excls = controlDict[pfx+'Exclude']
2420	for excl in excls:
2421	print >>pFile,' Excluded region: %8.2f%s to %8.2f%s'%(excl[0],units,excl[1],units)
2422	PrintSampleParms(Sample)
2423	PrintInstParms(Inst)
2424	PrintBackground(Background)
2425	elif 'HKLF' in histogram:
2426	Histogram = Histograms[histogram]
2427	hId = Histogram['hId']
2428	pfx = ':'+str(hId)+':'
2429	controlDict[pfx+'wtFactor'] = Histogram['wtFactor']
2430	Inst = Histogram['Instrument Parameters'][0]
2431	controlDict[pfx+'histType'] = Inst['Type'][0]
2432	if 'X' in Inst['Type'][0]:
2433	histDict[pfx+'Lam'] = Inst['Lam'][1]
2434	controlDict[pfx+'keV'] = 12.397639/histDict[pfx+'Lam']
2435	return histVary,histDict,controlDict
2436
2437	def SetHistogramData(parmDict,sigDict,Histograms,Print=True,pFile=None):
2438	'needs a doc string'
2439
2440	def SetBackgroundParms(pfx,Background,parmDict,sigDict):
2441	Back = Background[0]
2442	DebyePeaks = Background[1]
2443	lenBack = len(Back[3:])
2444	backSig = [0 for i in range(lenBack+3DebyePeaks['nDebye']+4DebyePeaks['nPeaks'])]
2445	for i in range(lenBack):
2446	Back[3+i] = parmDict[pfx+'Back;'+str(i)]
2447	if pfx+'Back;'+str(i) in sigDict:
2448	backSig[i] = sigDict[pfx+'Back;'+str(i)]
2449	if DebyePeaks['nDebye']:
2450	for i in range(DebyePeaks['nDebye']):
2451	names = [pfx+'DebyeA:'+str(i),pfx+'DebyeR:'+str(i),pfx+'DebyeU:'+str(i)]
2452	for j,name in enumerate(names):
2453	DebyePeaks['debyeTerms'][i][2*j] = parmDict[name]
2454	if name in sigDict:
2455	backSig[lenBack+3*i+j] = sigDict[name]
2456	if DebyePeaks['nPeaks']:
2457	for i in range(DebyePeaks['nPeaks']):
2458	names = [pfx+'BkPkpos;'+str(i),pfx+'BkPkint;'+str(i),
2459	pfx+'BkPksig;'+str(i),pfx+'BkPkgam;'+str(i)]
2460	for j,name in enumerate(names):
2461	DebyePeaks['peaksList'][i][2*j] = parmDict[name]
2462	if name in sigDict:
2463	backSig[lenBack+3DebyePeaks['nDebye']+4i+j] = sigDict[name]
2464	return backSig
2465
2466	def SetInstParms(pfx,Inst,parmDict,sigDict):
2467	instSig = {}
2468	insKeys = Inst.keys()
2469	insKeys.sort()
2470	for item in insKeys:
2471	insName = pfx+item
2472	Inst[item][1] = parmDict[insName]
2473	if insName in sigDict:
2474	instSig[item] = sigDict[insName]
2475	else:
2476	instSig[item] = 0
2477	return instSig
2478
2479	def SetSampleParms(pfx,Sample,parmDict,sigDict):
2480	if 'Bragg' in Sample['Type']: #Bragg-Brentano
2481	sampSig = [0 for i in range(5)]
2482	for i,item in enumerate(['Scale','Shift','Transparency','SurfRoughA','SurfRoughB']):
2483	Sample[item][0] = parmDict[pfx+item]
2484	if pfx+item in sigDict:
2485	sampSig[i] = sigDict[pfx+item]
2486	elif 'Debye' in Sample['Type']: #Debye-Scherrer
2487	sampSig = [0 for i in range(4)]
2488	for i,item in enumerate(['Scale','Absorption','DisplaceX','DisplaceY']):
2489	Sample[item][0] = parmDict[pfx+item]
2490	if pfx+item in sigDict:
2491	sampSig[i] = sigDict[pfx+item]
2492	return sampSig
2493
2494	def PrintBackgroundSig(Background,backSig):
2495	Back = Background[0]
2496	DebyePeaks = Background[1]
2497	lenBack = len(Back[3:])
2498	valstr = ' value : '
2499	sigstr = ' sig : '
2500	refine = False
2501	for i,back in enumerate(Back[3:]):
2502	valstr += '%10.4g'%(back)
2503	if Back[1]:
2504	refine = True
2505	sigstr += '%10.4g'%(backSig[i])
2506	else:
2507	sigstr += 10*' '
2508	if refine:
2509	print >>pFile,'\n Background function: ',Back[0]
2510	print >>pFile,valstr
2511	print >>pFile,sigstr
2512	if DebyePeaks['nDebye']:
2513	ifAny = False
2514	ptfmt = "%12.3f"
2515	names = ' names :'
2516	ptstr = ' values:'
2517	sigstr = ' esds :'
2518	for item in sigDict:
2519	if 'Debye' in item:
2520	ifAny = True
2521	names += '%12s'%(item)
2522	ptstr += ptfmt%(parmDict[item])
2523	sigstr += ptfmt%(sigDict[item])
2524	if ifAny:
2525	print >>pFile,'\n Debye diffuse scattering coefficients'
2526	print >>pFile,names
2527	print >>pFile,ptstr
2528	print >>pFile,sigstr
2529	if DebyePeaks['nPeaks']:
2530	ifAny = False
2531	ptfmt = "%14.3f"
2532	names = ' names :'
2533	ptstr = ' values:'
2534	sigstr = ' esds :'
2535	for item in sigDict:
2536	if 'BkPk' in item:
2537	ifAny = True
2538	names += '%14s'%(item)
2539	ptstr += ptfmt%(parmDict[item])
2540	sigstr += ptfmt%(sigDict[item])
2541	if ifAny:
2542	print >>pFile,'\n Single peak coefficients'
2543	print >>pFile,names
2544	print >>pFile,ptstr
2545	print >>pFile,sigstr
2546
2547	def PrintInstParmsSig(Inst,instSig):
2548	refine = False
2549	insKeys = instSig.keys()
2550	insKeys.sort()
2551	iBeg = 0
2552	Ok = True
2553	while Ok:
2554	ptlbls = ' names :'
2555	ptstr = ' value :'
2556	sigstr = ' sig :'
2557	iFin = min(iBeg+9,len(insKeys))
2558	for name in insKeys[iBeg:iFin]:
2559	if name not in ['Type','Lam1','Lam2','Azimuth','Source','fltPath']:
2560	ptlbls += '%12s' % (name)
2561	ptstr += '%12.6f' % (Inst[name][1])
2562	if instSig[name]:
2563	refine = True
2564	sigstr += '%12.6f' % (instSig[name])
2565	else:
2566	sigstr += 12*' '
2567	if refine:
2568	print >>pFile,'\n Instrument Parameters:'
2569	print >>pFile,ptlbls
2570	print >>pFile,ptstr
2571	print >>pFile,sigstr
2572	iBeg = iFin
2573	if iBeg == len(insKeys):
2574	Ok = False
2575
2576	def PrintSampleParmsSig(Sample,sampleSig):
2577	ptlbls = ' names :'
2578	ptstr = ' values:'
2579	sigstr = ' sig :'
2580	refine = False
2581	if 'Bragg' in Sample['Type']:
2582	for i,item in enumerate(['Scale','Shift','Transparency','SurfRoughA','SurfRoughB']):
2583	ptlbls += '%14s'%(item)
2584	ptstr += '%14.4f'%(Sample[item][0])
2585	if sampleSig[i]:
2586	refine = True
2587	sigstr += '%14.4f'%(sampleSig[i])
2588	else:
2589	sigstr += 14*' '
2590
2591	elif 'Debye' in Sample['Type']: #Debye-Scherrer
2592	for i,item in enumerate(['Scale','Absorption','DisplaceX','DisplaceY']):
2593	ptlbls += '%14s'%(item)
2594	ptstr += '%14.4f'%(Sample[item][0])
2595	if sampleSig[i]:
2596	refine = True
2597	sigstr += '%14.4f'%(sampleSig[i])
2598	else:
2599	sigstr += 14*' '
2600
2601	if refine:
2602	print >>pFile,'\n Sample Parameters:'
2603	print >>pFile,ptlbls
2604	print >>pFile,ptstr
2605	print >>pFile,sigstr
2606
2607	histoList = Histograms.keys()
2608	histoList.sort()
2609	for histogram in histoList:
2610	if 'PWDR' in histogram:
2611	Histogram = Histograms[histogram]
2612	hId = Histogram['hId']
2613	pfx = ':'+str(hId)+':'
2614	Background = Histogram['Background']
2615	backSig = SetBackgroundParms(pfx,Background,parmDict,sigDict)
2616
2617	Inst = Histogram['Instrument Parameters'][0]
2618	instSig = SetInstParms(pfx,Inst,parmDict,sigDict)
2619
2620	Sample = Histogram['Sample Parameters']
2621	sampSig = SetSampleParms(pfx,Sample,parmDict,sigDict)
2622
2623	print >>pFile,'\n Histogram: ',histogram,' histogram Id: ',hId
2624	print >>pFile,135*'-'
2625	print >>pFile,' PWDR histogram weight factor = '+'%.3f'%(Histogram['wtFactor'])
2626	print >>pFile,' Final refinement wR = %.2f%% on %d observations in this histogram'%(Histogram['Residuals']['wR'],Histogram['Residuals']['Nobs'])
2627	print >>pFile,' Other residuals: R = %.2f%%, Rb = %.2f%%, wRb = %.2f%% wRmin = %.2f%%'% \
2628	(Histogram['Residuals']['R'],Histogram['Residuals']['Rb'],Histogram['Residuals']['wRb'],Histogram['Residuals']['wRmin'])
2629	if Print:
2630	print >>pFile,' Instrument type: ',Sample['Type']
2631	PrintSampleParmsSig(Sample,sampSig)
2632	PrintInstParmsSig(Inst,instSig)
2633	PrintBackgroundSig(Background,backSig)
2634

Note: See TracBrowser for help on using the repository browser.

Download in other formats: