sans_gen.py @ d7b6913

ESS_GUIESS_GUI_DocsESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_iss959ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalccostrafo411magnetic_scattrelease-4.1.1release-4.1.2release-4.2.2release_4.0.1ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload

Last change on this file since d7b6913 was d7b6913, checked in by Jae Cho <jhjcho@…>, 11 years ago
sans gen utest
Property mode set to `100644`
File size: 36.5 KB

Line
1	"""
2	SANS generic computation and sld file readers
3	"""
4	from sans.models.BaseComponent import BaseComponent
5	import sans.models.sans_extension.sld2i as mod
6	from periodictable import formula
7	import numpy
8	import os
9
10	MFactor_AM = 2.853E-12
11	MFactor_mT = 2.3164E-9
12	METER2ANG = 1.0E+10
13
14	def mag2sld(mag, v_unit=None):
15	"""
16	Convert magnetization to magnatic SLD
17	sldm = Dm * mag where Dm = gamma * classical elec. radius/(2*Bohr magneton)
18	Dm ~ 2.853E-12 [A^(-2)] ==> Shouldn't be 2.90636E-12 [A^(-2)]???
19	"""
20	if v_unit == "A/m":
21	factor = MFactor_AM
22	elif v_unit == "mT":
23	factor = MFactor_mT
24	else:
25	raise ValueError, "Invalid valueunit"
26	sld_m = factor * mag
27	return sld_m
28
29	class GenSAS(BaseComponent):
30	"""
31	Generic SAS computation Model based on sld (n & m) arrays
32	"""
33	def __init__(self):
34	"""
35	Init
36
37	:Params sld_data: MagSLD object
38	"""
39	# Initialize BaseComponent
40	BaseComponent.__init__(self)
41	self.sld_data = None
42	self.data_pos_unit = None
43	self.data_x = None
44	self.data_y = None
45	self.data_z = None
46	self.data_sldn = None
47	self.data_mx = None
48	self.data_my = None
49	self.data_mz = None
50	self.volume = 1000 #[A^3]
51	## Name of the model
52	self.name = "GenSAS"
53	## Define parameters
54	self.params = {}
55	self.params['scale'] = 1.0
56	self.params['background'] = 0.0
57	self.params['Up_frac_in'] = 0.5
58	self.params['Up_frac_out'] = 0.5
59	self.params['Up_theta'] = 0.0
60	self.description='GenSAS'
61	## Parameter details [units, min, max]
62	self.details = {}
63	self.details['scale'] = ['', None, None]
64	self.details['background'] = ['[1/cm]', None, None]
65	self.details['Up_frac_in'] = ['[u/(u+d)]', None, None]
66	self.details['Up_frac_out'] = ['[u/(u+d)]', None, None]
67	self.details['Up_theta'] = ['[deg]', None, None]
68	# fixed parameters
69	self.fixed=[]
70
71	def set_pixel_volume(self, volume):
72	"""
73	Set the volume of a pixel in (A^3) unit
74	"""
75	self.volume = volume
76
77	def _gen(self, x, y, i):
78	"""
79	Evaluate the function
80	:Param x: array of x-values
81	:Param y: array of y-values
82	:Param i: array of initial i-value
83	:return: function value
84	"""
85	len_x = len(self.data_x)
86	len_q = len(x)
87	model = mod.new_GenI(len_x, self.data_x, self.data_y, self.data_z,
88	self.data_sldn, self.data_mx, self.data_my,
89	self.data_mz, self.params['Up_frac_in'],
90	self.params['Up_frac_out'],
91	self.params['Up_theta'])
92	mod.genicom(model, len_q, x, y, i)
93	return self.params['scale'] * i * self.volume \
94	+ self.params['background']
95
96	def set_sld_data(self, sld_data=None):
97	"""
98	Sets sld_data
99	"""
100	self.sld_data = sld_data
101	self.data_pos_unit = sld_data.pos_unit
102	self.data_x = sld_data.pos_x
103	self.data_y = sld_data.pos_y
104	self.data_z = sld_data.pos_z
105	self.data_sldn = sld_data.sld_n
106	self.data_mx = sld_data.sld_mx
107	self.data_my = sld_data.sld_my
108	self.data_mz = sld_data.sld_mz
109
110	def getProfile(self):
111	"""
112	Get SLD profile
113
114	: return: sld_data
115	"""
116	return self.sld_data
117
118	def run(self, x = 0.0):
119	"""
120	Evaluate the model
121	:param x: simple value
122	:return: (I value)
123	"""
124	if x.__class__.__name__ == 'list':
125	i_out = numpy.zeros_like(x[0])
126	y_in = numpy.zero_like(x[0])
127	# 1D I is found at y =0 in the 2D pattern
128	return self._gen(x[0],y_in, i_out )
129	else:
130	msg = "Q must be given as list of qx's and qy's"
131	raise ValueError, msg
132
133	def runXY(self, x = 0.0):
134	"""
135	Evaluate the model
136	:param x: simple value
137	:return: I value
138	:Use this runXY() for the computation
139	"""
140	if x.__class__.__name__ == 'list':
141	i_out = numpy.zeros_like(x[0])
142	import time
143	s=time.time()
144	out = self._gen(x[0], x[1], i_out)
145	print "i_out_time",time.time()-s
146	return out
147	else:
148	msg = "Q must be given as list of qx's and qy's"
149	raise ValueError, msg
150
151	def evalDistribution(self, qdist):
152	"""
153	Evaluate a distribution of q-values.
154	* For 1D, a numpy array is expected as input:
155	evalDistribution(q)
156	where q is a numpy array.
157	* For 2D, a list of numpy arrays are expected: [qx_prime,qy_prime],
158	where 1D arrays,
159	:param qdist: ndarray of scalar q-values or list [qx,qy]
160	where qx,qy are 1D ndarrays
161	"""
162	if qdist.__class__.__name__ == 'list':
163	out = self.runXY(qdist)
164	return out
165	else:
166	mesg = "evalDistribution is expecting an ndarray of "
167	mesg += "a list [qx,qy] where qx,qy are arrays."
168	raise RuntimeError, mesg
169
170	class OMF2SLD:
171	"""
172	Convert OMFData to MAgData
173	"""
174	def __init__(self):
175	"""
176	Init
177	"""
178	self.pos_x = None
179	self.pos_y = None
180	self.pos_z = None
181	self.mx = None
182	self.my = None
183	self.mz = None
184	self.sld_n = None
185	self.output = None
186	self.omfdata = None
187
188	def set_data(self, omfdata, shape='rectangular'):
189	"""
190	Set all data
191	"""
192	self.omfdata = omfdata
193	length = int(omfdata.xnodes * omfdata.ynodes * omfdata.znodes)
194	pos_x = numpy.arange(omfdata.xmin,
195	omfdata.xnodes*omfdata.xstepsize + omfdata.xmin,
196	omfdata.xstepsize)
197	pos_y = numpy.arange(omfdata.ymin,
198	omfdata.ynodes*omfdata.ystepsize + omfdata.ymin,
199	omfdata.ystepsize)
200	pos_z = numpy.arange(omfdata.zmin,
201	omfdata.znodes*omfdata.zstepsize + omfdata.zmin,
202	omfdata.zstepsize)
203	self.pos_x = numpy.tile(pos_x, int(omfdata.ynodes * omfdata.znodes))
204	self.pos_y = pos_y.repeat(int(omfdata.xnodes))
205	self.pos_y = numpy.tile(self.pos_y, int(omfdata.znodes))
206	self.pos_z = pos_z.repeat(int(omfdata.xnodes * omfdata.ynodes))
207	self.mx = omfdata.mx
208	self.my = omfdata.my
209	self.mz = omfdata.mz
210	self.sld_n = numpy.zeros(length)
211
212	if omfdata.mx == None:
213	self.mx = numpy.zeros(length)
214	if omfdata.my == None:
215	self.my = numpy.zeros(length)
216	if omfdata.mz == None:
217	self.mz = numpy.zeros(length)
218
219	self._check_data_length(length)
220	self.remove_null_points()
221	mask = numpy.ones(len(self.sld_n), dtype=bool)
222	if shape == 'ellipsoid':
223	try:
224	# Pixel (step) size included
225	x_r = (max(self.pos_x) - min(self.pos_x) + \
226	omfdata.xstepsize) / 2.0
227	y_r = (max(self.pos_y) - min(self.pos_y) + \
228	omfdata.ystepsize) / 2.0
229	z_r = (max(self.pos_z) - min(self.pos_z) + \
230	omfdata.zstepsize) / 2.0
231	x_dir2 = (self.pos_x / x_r) * (self.pos_x / x_r)
232	y_dir2 = (self.pos_y / y_r) * (self.pos_y / y_r)
233	z_dir2 = (self.pos_z / z_r) * (self.pos_z / z_r)
234	mask = (x_dir2 + y_dir2 + z_dir2) <= 1
235	except:
236	pass
237	self.output = MagSLD(self.pos_x[mask], self.pos_y[mask],
238	self.pos_z[mask], self.sld_n[mask],
239	self.mx[mask], self.my[mask], self.mz[mask])
240	self.output.set_pix_type('pixel')
241	self.output.set_pixel_symbols('pixel')
242
243	def get_omfdata(self):
244	"""
245	Return all data
246	"""
247	return self.omfdata
248
249	def get_output(self):
250	"""
251	Return output
252	"""
253	return self.output
254
255	def _check_data_length(self, length):
256	"""
257	Check if the data lengths are consistent
258	:Params length: data length
259	"""
260	msg = "Error: Inconsistent data length."
261	if len(self.pos_x) != length:
262	raise ValueError, msg
263	if len(self.pos_y) != length:
264	raise ValueError, msg
265	if len(self.pos_z) != length:
266	raise ValueError, msg
267	if len(self.mx) != length:
268	raise ValueError, msg
269	if len(self.my) != length:
270	raise ValueError, msg
271	if len(self.mz) != length:
272	raise ValueError, msg
273
274	def remove_null_points(self, remove=False):
275	"""
276	Removes any mx, my, and mz = 0 points
277	"""
278	if remove:
279	is_nonzero = (numpy.fabs(self.mx) + numpy.fabs(self.my) +
280	numpy.fabs(self.mz)).nonzero()
281	if len(is_nonzero[0]) > 0:
282	self.pos_x = self.pos_x[is_nonzero]
283	self.pos_y = self.pos_y[is_nonzero]
284	self.pos_z = self.pos_z[is_nonzero]
285	self.sld_n = self.sld_n[is_nonzero]
286	self.mx = self.mx[is_nonzero]
287	self.my = self.my[is_nonzero]
288	self.mz = self.mz[is_nonzero]
289	self.pos_x -= (min(self.pos_x) + max(self.pos_x)) / 2.0
290	self.pos_y -= (min(self.pos_y) + max(self.pos_y)) / 2.0
291	self.pos_z -= (min(self.pos_z) + max(self.pos_z)) / 2.0
292
293	def get_magsld(self):
294	"""
295	return MagSLD
296	"""
297	return self.output
298
299
300	class OMFReader:
301	"""
302	Class to load omf/ascii files (3 columns w/header).
303	"""
304	## File type
305	type_name = "OMF ASCII"
306
307	## Wildcards
308	type = ["OMF files (.OMF, .omf)\|*.omf"]
309	## List of allowed extensions
310	ext = ['.omf', '.OMF']
311
312	def read(self, path):
313	"""
314	Load data file
315	:param path: file path
316	:return: x, y, z, sld_n, sld_mx, sld_my, sld_mz
317	"""
318	data_conv_m = None
319	desc = ""
320	mx = numpy.zeros(0)
321	my = numpy.zeros(0)
322	mz = numpy.zeros(0)
323	try:
324	data_conv_val = None
325	data_conv_mesh = None
326	input_f = open(path, 'rb')
327	buff = input_f.read()
328	lines = buff.split('\n')
329	input_f.close()
330	output = OMFData()
331	valueunit = None
332	for line in lines:
333	toks = line.split()
334	# Read data
335	try:
336	_mx = float(toks[0])
337	_my = float(toks[1])
338	_mz = float(toks[2])
339	_mx = mag2sld(_mx, valueunit)
340	_my = mag2sld(_my, valueunit)
341	_mz = mag2sld(_mz, valueunit)
342	mx = numpy.append(mx, _mx)
343	my = numpy.append(my, _my)
344	mz = numpy.append(mz, _mz)
345	except:
346	# Skip non-data lines
347	pass
348	#Reading Header; Segment count ignored
349	s_line = line.split(":", 1)
350	if s_line[0].lower().count("oommf") > 0:
351	oommf = s_line[1].lstrip()
352	if s_line[0].lower().count("title") > 0:
353	title = s_line[1].lstrip()
354	if s_line[0].lower().count("desc") > 0:
355	desc += s_line[1].lstrip()
356	desc += '\n'
357	if s_line[0].lower().count("meshtype") > 0:
358	meshtype = s_line[1].lstrip()
359	if s_line[0].lower().count("meshunit") > 0:
360	meshunit = s_line[1].lstrip()
361	if meshunit.count("m") < 1:
362	msg = "Error: \n"
363	msg += "We accept only m as meshunit"
364	raise ValueError, msg
365	if s_line[0].lower().count("xbase") > 0:
366	xbase = s_line[1].lstrip()
367	if s_line[0].lower().count("ybase") > 0:
368	ybase = s_line[1].lstrip()
369	if s_line[0].lower().count("zbase") > 0:
370	zbase = s_line[1].lstrip()
371	if s_line[0].lower().count("xstepsize") > 0:
372	xstepsize = s_line[1].lstrip()
373	if s_line[0].lower().count("ystepsize") > 0:
374	ystepsize = s_line[1].lstrip()
375	if s_line[0].lower().count("zstepsize") > 0:
376	zstepsize = s_line[1].lstrip()
377	if s_line[0].lower().count("xnodes") > 0:
378	xnodes = s_line[1].lstrip()
379	if s_line[0].lower().count("ynodes") > 0:
380	ynodes = s_line[1].lstrip()
381	if s_line[0].lower().count("znodes") > 0:
382	znodes = s_line[1].lstrip()
383	if s_line[0].lower().count("xmin") > 0:
384	xmin = s_line[1].lstrip()
385	if s_line[0].lower().count("ymin") > 0:
386	ymin = s_line[1].lstrip()
387	if s_line[0].lower().count("zmin") > 0:
388	zmin = s_line[1].lstrip()
389	if s_line[0].lower().count("xmax") > 0:
390	xmax = s_line[1].lstrip()
391	if s_line[0].lower().count("ymax") > 0:
392	ymax = s_line[1].lstrip()
393	if s_line[0].lower().count("zmax") > 0:
394	zmax = s_line[1].lstrip()
395	if s_line[0].lower().count("valueunit") > 0:
396	valueunit = s_line[1].lstrip().rstrip()
397	if s_line[0].lower().count("valuemultiplier") > 0:
398	valuemultiplier = s_line[1].lstrip()
399	if s_line[0].lower().count("valuerangeminmag") > 0:
400	valuerangeminmag = s_line[1].lstrip()
401	if s_line[0].lower().count("valuerangemaxmag") > 0:
402	valuerangemaxmag = s_line[1].lstrip()
403	if s_line[0].lower().count("end") > 0:
404	#output.set_sldms(mx, my, mz)
405	output.filename = os.path.basename(path)
406	output.oommf = oommf
407	output.title = title
408	output.desc = desc
409	output.meshtype = meshtype
410	output.xbase = float(xbase) * METER2ANG
411	output.ybase = float(ybase) * METER2ANG
412	output.zbase = float(zbase) * METER2ANG
413	output.xstepsize = float(xstepsize) * METER2ANG
414	output.ystepsize = float(ystepsize) * METER2ANG
415	output.zstepsize = float(zstepsize) * METER2ANG
416	output.xnodes = float(xnodes)
417	output.ynodes = float(ynodes)
418	output.znodes = float(znodes)
419	output.xmin = float(xmin) * METER2ANG
420	output.ymin = float(ymin) * METER2ANG
421	output.zmin = float(zmin) * METER2ANG
422	output.xmax = float(xmax) * METER2ANG
423	output.ymax = float(ymax) * METER2ANG
424	output.zmax = float(zmax) * METER2ANG
425	output.valuemultiplier = valuemultiplier
426	output.valuerangeminmag = mag2sld(float(valuerangeminmag),\
427	valueunit)
428	output.valuerangemaxmag = mag2sld(float(valuerangemaxmag),\
429	valueunit)
430	output.set_m(mx, my, mz)
431	return output
432	except:
433	msg = "%s is not supported: \n" % path
434	msg += "We accept only Text format OMF file."
435	raise RuntimeError, msg
436
437	class PDBReader:
438	"""
439	PDB reader class: limited for reading the lines starting with 'ATOM'
440	"""
441	type_name = "PDB"
442	## Wildcards
443	type = ["pdb files (.PDB, .pdb)\|*.pdb"]
444	## List of allowed extensions
445	ext = ['.pdb', '.PDB']
446
447	def read(self, path):
448	"""
449	Load data file
450
451	:param path: file path
452
453	:return: MagSLD
454
455	:raise RuntimeError: when the file can't be opened
456	"""
457	pos_x = numpy.zeros(0)
458	pos_y = numpy.zeros(0)
459	pos_z = numpy.zeros(0)
460	sld_n = numpy.zeros(0)
461	sld_mx = numpy.zeros(0)
462	sld_my = numpy.zeros(0)
463	sld_mz = numpy.zeros(0)
464	pix_symbol = numpy.zeros(0)
465	try:
466	input_f = open(path, 'rb')
467	buff = input_f.read()
468	lines = buff.split('\n')
469	input_f.close()
470	for line in lines:
471	try:
472	#toks = line.split()
473	# check if line starts with "ATOM"
474	if line[0:6].strip().count('ATM') > 0 or \
475	line[0:6].strip() == 'ATOM':
476	# define fields of interest
477	atom_id = line[6:11].strip()
478	atom_name = line[12:16].strip()
479	try:
480	float(line[12])
481	atom_name = atom_name[1].upper()
482	except:
483	if len(atom_name) == 4:
484	atom_name = atom_name[0].upper()
485	elif line[12] != ' ':
486	atom_name = atom_name[0].upper() + \
487	atom_name[1].lower()
488	else:
489	atom_name = atom_name[0].upper()
490	#res_name = line[17:20].strip()
491	#chain_name = line[21:22].strip()
492	_pos_x = float(line[30:38].strip())
493	_pos_y = float(line[38:46].strip())
494	_pos_z = float(line[46:54].strip())
495
496	pos_x = numpy.append(pos_x, _pos_x)
497	pos_y = numpy.append(pos_y, _pos_y)
498	pos_z = numpy.append(pos_z, _pos_z)
499	try:
500	#from periodictable import nsf
501	#val = nsf.neutron_sld(atom_name)
502	val = formula(atom_name).atoms.keys()[0].neutron.b_c
503	# sld in Ang unit (from fm)
504	val = float(val)
505	val *= 1.0e-5
506	sld_n = numpy.append(sld_n, val)
507	except:
508	print "Error: set the sld of %s to zero"% atom_name
509	sld_n = numpy.append(sld_n, 0.0)
510	sld_mx = numpy.append(sld_mx, 0)
511	sld_my = numpy.append(sld_my, 0)
512	sld_mz = numpy.append(sld_mz, 0)
513	pix_symbol = numpy.append(pix_symbol, atom_name)
514	except:
515	pass
516	#print pix_symbol, pos_x
517	pos_x -= (min(pos_x) + max(pos_x)) / 2.0
518	pos_y -= (min(pos_y) + max(pos_y)) / 2.0
519	pos_z -= (min(pos_z) + max(pos_z)) / 2.0
520
521	output = MagSLD(pos_x, pos_y, pos_z, sld_n, sld_mx, sld_my, sld_mz)
522	#output.set_sldms(0, 0, 0)
523	output.filename = os.path.basename(path)
524	output.set_pix_type('atom')
525	output.set_pixel_symbols(pix_symbol)
526	output.set_nodes()
527	output.sld_unit = 'A'
528	return output
529	except:
530	RuntimeError, "%s is not a sld file" % path
531
532	def write(self, path, data):
533	"""
534	Write
535	"""
536	#Not implemented
537	print "Not implemented... "
538
539	class SLDReader:
540	"""
541	Class to load ascii files (7 columns).
542	"""
543	## File type
544	type_name = "SLD ASCII"
545
546	## Wildcards
547	type = ["sld files (.SLD, .sld)\|*.sld",
548	"txt files (.TXT, .txt)\|*.txt"]
549	## List of allowed extensions
550	ext = ['.txt', '.TXT', '.sld', 'SLD', '.*']
551
552	def read(self, path):
553	"""
554	Load data file
555
556	:param path: file path
557
558	:return MagSLD: x, y, z, sld_n, sld_mx, sld_my, sld_mz
559
560	:raise RuntimeError: when the file can't be opened
561	:raise ValueError: when the length of the data vectors are inconsistent
562	"""
563	pos_x = numpy.zeros(0)
564	pos_y = numpy.zeros(0)
565	pos_z = numpy.zeros(0)
566	sld_n = numpy.zeros(0)
567	sld_mx = numpy.zeros(0)
568	sld_my = numpy.zeros(0)
569	sld_mz = numpy.zeros(0)
570	try:
571	input_f = open(path, 'rb')
572	buff = input_f.read()
573	lines = buff.split('\n')
574	input_f.close()
575	for line in lines:
576	toks = line.split()
577	try:
578	_pos_x = float(toks[0])
579	_pos_y = float(toks[1])
580	_pos_z = float(toks[2])
581	_sld_n = float(toks[3])
582	_sld_mx = float(toks[4])
583	_sld_my = float(toks[5])
584	_sld_mz = float(toks[6])
585
586	pos_x = numpy.append(pos_x, _pos_x)
587	pos_y = numpy.append(pos_y, _pos_y)
588	pos_z = numpy.append(pos_z, _pos_z)
589	sld_n = numpy.append(sld_n, _sld_n)
590	sld_mx = numpy.append(sld_mx, _sld_mx)
591	sld_my = numpy.append(sld_my, _sld_my)
592	sld_mz = numpy.append(sld_mz, _sld_mz)
593	except:
594	# Skip non-data lines
595	pass
596	pos_x -= (min(pos_x) + max(pos_x)) / 2.0
597	pos_y -= (min(pos_y) + max(pos_y)) / 2.0
598	pos_z -= (min(pos_z) + max(pos_z)) / 2.0
599
600	output = MagSLD(pos_x, pos_y, pos_z, sld_n, sld_mx, sld_my, sld_mz)
601	output.filename = os.path.basename(path)
602	output.set_pix_type('pixel')
603	output.set_pixel_symbols('pixel')
604	return output
605	except:
606	RuntimeError, "%s is not a sld file" % path
607
608	def write(self, path, data):
609	"""
610	Write sld file
611
612	:Param path: file path
613	:Param data: MagSLD data object
614	"""
615	if path == None:
616	raise ValueError, "Missing the file path."
617	if data == None:
618	raise ValueError, "Missing the data to save."
619
620	x_val = data.pos_x
621	y_val = data.pos_y
622	z_val = data.pos_z
623
624	length = len(x_val)
625
626	sld_n = data.sld_n
627	if sld_n == None:
628	sld_n = numpy.zerros(length)
629
630	sld_mx = data.sld_mx
631	if sld_mx == None:
632	sld_mx = numpy.zerros(length)
633	sld_my = numpy.zerros(length)
634	sld_mz = numpy.zerros(length)
635	else:
636	sld_my = data.sld_my
637	sld_mz = data.sld_mz
638
639	out = open(path, 'w')
640	# First Line: Column names
641	out.write("X Y Z SLDN SLDMx SLDMy SLDMz")
642	for ind in range(length):
643	out.write("\n%g %g %g %g %g %g %g" % (x_val[ind], y_val[ind],
644	z_val[ind], sld_n[ind],
645	sld_mx[ind], sld_my[ind],
646	sld_mz[ind]))
647	out.close()
648
649
650	class OMFData:
651	"""
652	OMF Data.
653	"""
654	_meshunit = "A"
655	_valueunit = "A^(-2)"
656	def __init__(self):
657	"""
658	Init for mag SLD
659	"""
660	self.filename = 'default'
661	self.oommf = ''
662	self.title = ''
663	self.desc = ''
664	self.meshtype = ''
665	self.meshunit = self._meshunit
666	self.valueunit = self._valueunit
667	self.xbase = 0.0
668	self.ybase = 0.0
669	self.zbase = 0.0
670	self.xstepsize = 6.0
671	self.ystepsize = 6.0
672	self.zstepsize = 6.0
673	self.xnodes = 10.0
674	self.ynodes = 10.0
675	self.znodes = 10.0
676	self.xmin = 0.0
677	self.ymin = 0.0
678	self.zmin = 0.0
679	self.xmax = 60.0
680	self.ymax = 60.0
681	self.zmax = 60.0
682	self.mx = None
683	self.my = None
684	self.mz = None
685	self.valuemultiplier = 1.
686	self.valuerangeminmag = 0
687	self.valuerangemaxmag = 0
688
689	def __str__(self):
690	"""
691	doc strings
692	"""
693	_str = "Type: %s\n" % self.__class__.__name__
694	_str += "File: %s\n" % self.filename
695	_str += "OOMMF: %s\n" % self.oommf
696	_str += "Title: %s\n" % self.title
697	_str += "Desc: %s\n" % self.desc
698	_str += "meshtype: %s\n" % self.meshtype
699	_str += "meshunit: %s\n" % str(self.meshunit)
700	_str += "xbase: %s [%s]\n" % (str(self.xbase), self.meshunit)
701	_str += "ybase: %s [%s]\n" % (str(self.ybase), self.meshunit)
702	_str += "zbase: %s [%s]\n" % (str(self.zbase), self.meshunit)
703	_str += "xstepsize: %s [%s]\n" % (str(self.xstepsize),
704	self.meshunit)
705	_str += "ystepsize: %s [%s]\n" % (str(self.ystepsize),
706	self.meshunit)
707	_str += "zstepsize: %s [%s]\n" % (str(self.zstepsize),
708	self.meshunit)
709	_str += "xnodes: %s\n" % str(self.xnodes)
710	_str += "ynodes: %s\n" % str(self.ynodes)
711	_str += "znodes: %s\n" % str(self.znodes)
712	_str += "xmin: %s [%s]\n" % (str(self.xmin), self.meshunit)
713	_str += "ymin: %s [%s]\n" % (str(self.ymin), self.meshunit)
714	_str += "zmin: %s [%s]\n" % (str(self.zmin), self.meshunit)
715	_str += "xmax: %s [%s]\n" % (str(self.xmax), self.meshunit)
716	_str += "ymax: %s [%s]\n" % (str(self.ymax), self.meshunit)
717	_str += "zmax: %s [%s]\n" % (str(self.zmax), self.meshunit)
718	_str += "valueunit: %s\n" % self.valueunit
719	_str += "valuemultiplier: %s\n" % str(self.valuemultiplier)
720	_str += "ValueRangeMinMag:%s [%s]\n" % (str(self.valuerangeminmag),
721	self.valueunit)
722	_str += "ValueRangeMaxMag:%s [%s]\n" % (str(self.valuerangemaxmag),
723	self.valueunit)
724	return _str
725
726	def set_m(self, mx, my, mz):
727	"""
728	Set the Mx, My, Mz values
729	"""
730	self.mx = mx
731	self.my = my
732	self.mz = mz
733
734	class MagSLD:
735	"""
736	Magnetic SLD.
737	"""
738	pos_x = None
739	pos_y = None
740	pos_z = None
741	sld_n = None
742	sld_mx = None
743	sld_my = None
744	sld_mz = None
745	# Units
746	_pos_unit = 'A'
747	_sld_unit = '1/A^(2)'
748	_pix_type = 'pixel'
749
750	def __init__(self, pos_x, pos_y, pos_z, sld_n=None,
751	sld_mx=None, sld_my=None, sld_mz=None):
752	"""
753	Init for mag SLD
754	:params : All should be numpy 1D array
755	"""
756	self.is_data = True
757	self.filename = ''
758	self.xstepsize = 6.0
759	self.ystepsize = 6.0
760	self.zstepsize = 6.0
761	self.xnodes = 10.0
762	self.ynodes = 10.0
763	self.znodes = 10.0
764	self.has_stepsize = False
765	self.pos_unit = self._pos_unit
766	self.sld_unit = self._sld_unit
767	self.pix_type = 'pixel'
768	self.pos_x = pos_x
769	self.pos_y = pos_y
770	self.pos_z = pos_z
771	self.sld_n = sld_n
772	self.sld_mx = sld_mx
773	self.sld_my = sld_my
774	self.sld_mz = sld_mz
775	self.sld_m = None
776	self.sld_phi = None
777	self.sld_theta = None
778	self.pix_symbol = None
779	if sld_mx != None and sld_my != None and sld_mz != None:
780	self.set_sldms(sld_mx, sld_my, sld_mz)
781	self.set_nodes()
782
783	def __str__(self):
784	"""
785	doc strings
786	"""
787	_str = "Type: %s\n" % self.__class__.__name__
788	_str += "File: %s\n" % self.filename
789	_str += "Axis_unit: %s\n" % self.pos_unit
790	_str += "SLD_unit: %s\n" % self.sld_unit
791	return _str
792
793	def set_pix_type(self, pix_type):
794	"""
795	Set pixel type
796	:Param pix_type: string, 'pixel' or 'atom'
797	"""
798	self.pix_type = pix_type
799
800	def set_sldn(self, sld_n):
801	"""
802	Sets neutron SLD
803	"""
804	if sld_n.__class__.__name__ == 'float':
805	if self.is_data:
806	# For data, put the value to only the pixels w non-zero M
807	is_nonzero = (numpy.fabs(self.sld_mx) +
808	numpy.fabs(self.sld_my) +
809	numpy.fabs(self.sld_mz)).nonzero()
810	self.sld_n = numpy.zeros(len(self.pos_x))
811	if len(self.sld_n[is_nonzero]) > 0:
812	self.sld_n[is_nonzero] = sld_n
813	else:
814	self.sld_n.fill(sld_n)
815	else:
816	# For non-data, put the value to all the pixels
817	self.sld_n = numpy.ones(len(self.pos_x)) * sld_n
818	else:
819	self.sld_n = sld_n
820
821	def set_sldms(self, sld_mx, sld_my, sld_mz):
822	"""
823	Sets (\|m\|, m_theta, m_phi)
824	"""
825	if sld_mx.__class__.__name__ == 'float':
826	self.sld_mx = numpy.ones(len(self.pos_x)) * sld_mx
827	else:
828	self.sld_mx = sld_mx
829	if sld_my.__class__.__name__ == 'float':
830	self.sld_my = numpy.ones(len(self.pos_x)) * sld_my
831	else:
832	self.sld_my = sld_my
833	if sld_mz.__class__.__name__ == 'float':
834	self.sld_mz = numpy.ones(len(self.pos_x)) * sld_mz
835	else:
836	self.sld_mz = sld_mz
837
838	sld_m = numpy.sqrt(sld_mx * sld_mx + sld_my * sld_my + \
839	sld_mz * sld_mz)
840	self.sld_m = sld_m
841
842	def set_pixel_symbols(self, symbol='pixel'):
843	"""
844	Set pixel
845	:Params pixel: str; pixel or atomic symbol, or array of strings
846	"""
847	if self.sld_n == None:
848	return
849	if symbol.__class__.__name__ == 'str':
850	self.pix_symbol = numpy.repeat(symbol, len(self.sld_n))
851	else:
852	self.pix_symbol = symbol
853
854	def get_sldn(self):
855	"""
856	Returns nuclear sld
857	"""
858	return self.sld_n
859
860	def get_sld_mxyz(self):
861	"""
862	Returns (sld_mx, sldmy, sld_mz)
863	"""
864	return (self.sld_mx, self.sldmy, self.sld_mz)
865
866	def get_sld_m(self):
867	"""
868	Returns (sldm, sld_theta, sld_phi)
869	"""
870	return (self.sldm, self.sld_theta, self.sld_phi)
871
872	def set_nodes(self):
873	"""
874	Set xnodes, ynodes, and znodes
875	"""
876	self.set_stepsize()
877	if self.pix_type == 'pixel':
878	try:
879	xdist = (max(self.pos_x) - min(self.pos_x)) / self.xstepsize
880	ydist = (max(self.pos_y) - min(self.pos_y)) / self.ystepsize
881	zdist = (max(self.pos_z) - min(self.pos_z)) / self.zstepsize
882	self.xnodes = int(xdist) + 1
883	self.ynodes = int(ydist) + 1
884	self.znodes = int(zdist) + 1
885	except:
886	self.xnodes = None
887	self.ynodes = None
888	self.znodes = None
889	else:
890	self.xnodes = None
891	self.ynodes = None
892	self.znodes = None
893
894	def set_stepsize(self):
895	"""
896	Set xtepsize, ystepsize, and zstepsize
897	"""
898	if self.pix_type == 'pixel':
899	try:
900	xpos_pre = self.pos_x[0]
901	ypos_pre = self.pos_y[0]
902	zpos_pre = self.pos_z[0]
903	for x_pos in self.pos_x:
904	if xpos_pre != x_pos:
905	self.xstepsize = numpy.fabs(x_pos - xpos_pre)
906	break
907	for y_pos in self.pos_y:
908	if ypos_pre != y_pos:
909	self.ystepsize = numpy.fabs(y_pos - ypos_pre)
910	break
911	for z_pos in self.pos_z:
912	if zpos_pre != z_pos:
913	self.zstepsize = numpy.fabs(z_pos - zpos_pre)
914	break
915	self.has_stepsize = True
916	except:
917	self.xstepsize = None
918	self.ystepsize = None
919	self.zstepsize = None
920	self.has_stepsize = False
921	else:
922	self.xstepsize = None
923	self.ystepsize = None
924	self.zstepsize = None
925	self.has_stepsize = True
926	return self.xstepsize, self.ystepsize, self.zstepsize
927
928
929	def test_load():
930	from danse.common.plottools.arrow3d import Arrow3D
931	import os
932	dir = os.path.abspath(os.path.curdir)
933	print dir
934	for i in range(6):
935	dir, _ = os.path.split(dir)
936	tfile = os.path.join(dir, "test", "CoreXY_ShellZ.txt")
937	ofile = os.path.join(dir, "test", "A_Raw_Example-1.omf")
938	if os.path.isfile(tfile):
939	tfpath = tfile
940	ofpath = ofile
941	break
942	reader = SLDReader()
943	oreader = OMFReader()
944	output = reader.read(tfpath)
945	ooutput = oreader.read(ofpath)
946	foutput = OMF2SLD()
947	foutput.set_data(ooutput)
948
949	import matplotlib.pyplot as plt
950	from mpl_toolkits.mplot3d import Axes3D
951	fig = plt.figure()
952	ax = fig.gca(projection='3d')
953	ax.plot(output.pos_x, output.pos_y, output.pos_z, '.', c="g",
954	alpha = 0.7, markeredgecolor='gray',rasterized=True)
955	gap = 7
956	max_mx = max(output.sld_mx)
957	max_my = max(output.sld_my)
958	max_mz = max(output.sld_mz)
959	max_m = max(max_mx,max_my,max_mz)
960	x2 = output.pos_x+output.sld_mx/max_m * gap
961	y2 = output.pos_y+output.sld_my/max_m * gap
962	z2 = output.pos_z+output.sld_mz/max_m * gap
963	x_arrow = numpy.column_stack((output.pos_x,x2))
964	y_arrow = numpy.column_stack((output.pos_y,y2))
965	z_arrow = numpy.column_stack((output.pos_z,z2))
966	unit_x2 = output.sld_mx / max_m
967	unit_y2 = output.sld_my / max_m
968	unit_z2 = output.sld_mz / max_m
969	color_x = numpy.fabs(unit_x2 * 0.8)
970	color_y = numpy.fabs(unit_y2 * 0.8)
971	color_z = numpy.fabs(unit_z2 * 0.8)
972	colors = numpy.column_stack((color_x, color_y, color_z))
973	a = Arrow3D(None, x_arrow,y_arrow,z_arrow, colors, mutation_scale=10, lw=1,
974	arrowstyle="->", color="y", alpha = 0.5)
975
976	ax.add_artist(a)
977	plt.show()
978
979	def test():
980	import os
981	dir = os.path.abspath(os.path.curdir)
982	for i in range(3):
983	dir, _ = os.path.split(dir)
984	#tfile = os.path.join(dir, "test", "C_Example_Converted.txt")
985	ofile = os.path.join(dir, "test", "A_Raw_Example-1.omf")
986	if os.path.isfile(ofile):
987	#tfpath = tfile
988	ofpath = ofile
989	break
990	#reader = SLDReader()
991	oreader = OMFReader()
992	#output = reader.read(tfpath)
993	ooutput = oreader.read(ofpath)
994	foutput = OMF2SLD()
995	foutput.set_data(ooutput)
996	writer = SLDReader()
997	writer.write(os.path.join(os.path.dirname(ofpath), "out.txt"),
998	foutput.output)
999	model = GenSAS()
1000	model.set_sld_data(foutput.output)
1001	x = numpy.arange(1000)/10000. + 1e-5
1002	y = numpy.arange(1000)/10000. + 1e-5
1003	z = numpy.arange(1000)/10000. + 1e-5
1004	i = numpy.zeros(1000)
1005	out = model.runXY([x,y,i])
1006
1007	if __name__ == "__main__":
1008	#test()
1009	test_load()

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

SasView

source: sasview/sanscalculator/src/sans/calculator/sans_gen.py @ d7b6913

Download in other formats: