Context Navigation

source: trunk/GSASIIstrIO.py @ 1174

Last change on this file since 1174 was 1160, checked in by toby, 12 years ago
finish ISODISPLACE fixes; improve show var window; improve help window; add refine checkbox for newvars in constraints display
Property svn:keywords set to `Date Author Revision URL Id`
File size: 113.2 KB

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

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

Download in other formats: