sans_gen.py @ dbc5dc9b

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

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

SasView

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

Download in other formats: