inversion_state.py @ 3a08f369

magnetic_scattrelease-4.2.2ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload

Last change on this file since 3a08f369 was a0e6b1b, checked in by lewis, 7 years ago
Esnure background value is saved/loaded when analysis is saved/loaded
Property mode set to `100644`
File size: 18.4 KB

Rev	Line
[959eb01]	1	"""
	2	Handling of P(r) inversion states
	3	"""
	4	################################################################################
	5	#This software was developed by the University of Tennessee as part of the
	6	#Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
	7	#project funded by the US National Science Foundation.
	8	#
	9	#See the license text in license.txt
	10	#
	11	#copyright 2009, University of Tennessee
	12	################################################################################
	13
	14
	15	import time
	16	import os
	17	import sys
	18	import logging
	19	from lxml import etree
	20	from sas.sasgui.guiframe.dataFitting import Data1D
	21	from sas.sascalc.dataloader.readers.cansas_reader import Reader as CansasReader
	22	from sas.sascalc.dataloader.readers.cansas_reader import get_content
	23
	24	logger = logging.getLogger(__name__)
	25
	26	PRNODE_NAME = 'pr_inversion'
	27	CANSAS_NS = "cansas1d/1.0"
	28
	29	# Translation of names between stored and object data
	30	## List of P(r) inversion inputs
	31	in_list = [["nterms", "nfunc"],
	32	["d_max", "d_max"],
	33	["alpha", "alpha"],
	34	["slit_width", "width"],
	35	["slit_height", "height"],
	36	["qmin", "qmin"],
	37	["qmax", "qmax"],
[a0e6b1b]	38	["estimate_bck", "estimate_bck"],
	39	["bck_value", "bck_value"]]
[959eb01]	40
	41	## List of P(r) inversion outputs
	42	out_list = [["elapsed", "elapsed"],
	43	["rg", "rg"],
	44	["iq0", "iq0"],
	45	["bck", "bck"],
	46	["chi2", "chi2"],
	47	["osc", "osc"],
	48	["pos", "pos"],
	49	["pos_err", "pos_err"],
	50	["alpha_estimate", "alpha_estimate"],
	51	["nterms_estimate", "nterms_estimate"]]
	52
	53	class InversionState(object):
	54	"""
	55	Class to hold the state information of the InversionControl panel.
	56	"""
	57	def __init__(self):
	58	"""
	59	Default values
	60	"""
	61	# Input
	62	self.file = None
	63	self.estimate_bck = False
	64	self.timestamp = time.time()
[a0e6b1b]	65	self.bck_value = 0.0
[959eb01]	66
	67	# Inversion parameters
	68	self.nfunc = None
	69	self.d_max = None
	70	self.alpha = None
	71
	72	# Slit parameters
	73	self.height = None
	74	self.width = None
	75
	76	# Q range
	77	self.qmin = None
	78	self.qmax = None
	79
	80	# Outputs
	81	self.elapsed = None
	82	self.rg = None
	83	self.iq0 = None
	84	self.bck = None
	85	self.chi2 = None
	86	self.osc = None
	87	self.pos = None
	88	self.pos_err = None
	89
	90	# Estimates
	91	self.alpha_estimate = None
	92	self.nterms_estimate = None
	93
	94	# Data
	95	self.q = None
	96	self.iq_obs = None
	97	self.iq_calc = None
	98
	99	# Coefficients
	100	self.coefficients = None
	101	self.covariance = None
	102
	103	def __str__(self):
	104	"""
	105	Pretty print
	106
	107	:return: string representing the state
	108
	109	"""
	110	state = "File: %s\n" % self.file
	111	state += "Timestamp: %s\n" % self.timestamp
	112	state += "Estimate bck: %s\n" % str(self.estimate_bck)
[a0e6b1b]	113	state += "Bck Value: %s\n" % str(self.bck_value)
[959eb01]	114	state += "No. terms: %s\n" % str(self.nfunc)
	115	state += "D_max: %s\n" % str(self.d_max)
	116	state += "Alpha: %s\n" % str(self.alpha)
	117
	118	state += "Slit height: %s\n" % str(self.height)
	119	state += "Slit width: %s\n" % str(self.width)
	120
	121	state += "Qmin: %s\n" % str(self.qmin)
	122	state += "Qmax: %s\n" % str(self.qmax)
	123
	124	state += "\nEstimates:\n"
	125	state += " Alpha: %s\n" % str(self.alpha_estimate)
	126	state += " Nterms: %s\n" % str(self.nterms_estimate)
	127
	128	state += "\nOutputs:\n"
	129	state += " Elapsed: %s\n" % str(self.elapsed)
	130	state += " Rg: %s\n" % str(self.rg)
	131	state += " I(q=0): %s\n" % str(self.iq0)
	132	state += " Bck: %s\n" % str(self.bck)
	133	state += " Chi^2: %s\n" % str(self.chi2)
	134	state += " Oscillation:%s\n" % str(self.osc)
	135	state += " Positive: %s\n" % str(self.pos)
	136	state += " 1-sigma pos:%s\n" % str(self.pos_err)
	137
	138	return state
	139
	140	def toXML(self, file="pr_state.prv", doc=None, entry_node=None):
	141	"""
	142	Writes the state of the InversionControl panel to file, as XML.
	143
	144	Compatible with standalone writing, or appending to an
	145	already existing XML document. In that case, the XML document
	146	is required. An optional entry node in the XML document
	147	may also be given.
	148
	149	:param file: file to write to
	150	:param doc: XML document object [optional]
	151	:param entry_node: XML node within the XML document at which
	152	we will append the data [optional]
	153
	154	"""
	155	#TODO: Get this to work
	156	from xml.dom.minidom import getDOMImplementation
	157
	158	# Check whether we have to write a standalone XML file
	159	if doc is None:
	160	impl = getDOMImplementation()
	161
	162	doc_type = impl.createDocumentType(PRNODE_NAME, "1.0", "1.0")
	163
	164	newdoc = impl.createDocument(None, PRNODE_NAME, doc_type)
	165	top_element = newdoc.documentElement
	166	else:
	167	# We are appending to an existing document
	168	newdoc = doc
	169	top_element = newdoc.createElement(PRNODE_NAME)
	170	if entry_node is None:
	171	newdoc.documentElement.appendChild(top_element)
	172	else:
	173	entry_node.appendChild(top_element)
	174
	175	attr = newdoc.createAttribute("version")
	176	attr.nodeValue = '1.0'
	177	top_element.setAttributeNode(attr)
	178
	179	# File name
	180	element = newdoc.createElement("filename")
	181	if self.file is not None:
	182	element.appendChild(newdoc.createTextNode(str(self.file)))
	183	else:
	184	element.appendChild(newdoc.createTextNode(str(file)))
	185	top_element.appendChild(element)
	186
	187	element = newdoc.createElement("timestamp")
	188	element.appendChild(newdoc.createTextNode(time.ctime(self.timestamp)))
	189	attr = newdoc.createAttribute("epoch")
	190	attr.nodeValue = str(self.timestamp)
	191	element.setAttributeNode(attr)
	192	top_element.appendChild(element)
	193
	194	# Inputs
	195	inputs = newdoc.createElement("inputs")
	196	top_element.appendChild(inputs)
	197
	198	for item in in_list:
	199	element = newdoc.createElement(item[0])
	200	element.appendChild(newdoc.createTextNode(str(getattr(self, item[1]))))
	201	inputs.appendChild(element)
	202
	203	# Outputs
	204	outputs = newdoc.createElement("outputs")
	205	top_element.appendChild(outputs)
	206
	207	for item in out_list:
	208	element = newdoc.createElement(item[0])
	209	element.appendChild(newdoc.createTextNode(str(getattr(self, item[1]))))
	210	outputs.appendChild(element)
	211
	212	# Save output coefficients and its covariance matrix
	213	element = newdoc.createElement("coefficients")
	214	element.appendChild(newdoc.createTextNode(str(self.coefficients)))
	215	outputs.appendChild(element)
	216	element = newdoc.createElement("covariance")
	217	element.appendChild(newdoc.createTextNode(str(self.covariance)))
	218	outputs.appendChild(element)
	219
	220	# Save the file
	221	if doc is None:
	222	fd = open(file, 'w')
	223	fd.write(newdoc.toprettyxml())
	224	fd.close()
	225	return None
	226	else:
	227	return newdoc
	228
	229	def fromXML(self, file=None, node=None):
	230	"""
	231	Load a P(r) inversion state from a file
	232
	233	:param file: .prv file
	234	:param node: node of a XML document to read from
	235
	236	"""
	237	if file is not None:
	238	msg = "InversionState no longer supports non-CanSAS"
	239	msg += " format for P(r) files"
	240	raise RuntimeError, msg
	241
	242	if node.get('version') and node.get('version') == '1.0':
	243
	244	# Get file name
	245	entry = get_content('ns:filename', node)
	246	if entry is not None:
	247	self.file = entry.text.strip()
	248
	249	# Get time stamp
	250	entry = get_content('ns:timestamp', node)
	251	if entry is not None and entry.get('epoch'):
	252	try:
	253	self.timestamp = float(entry.get('epoch'))
	254	except:
	255	msg = "InversionState.fromXML: Could not read "
	256	msg += "timestamp\n %s" % sys.exc_value
	257	logger.error(msg)
	258
	259	# Parse inversion inputs
	260	entry = get_content('ns:inputs', node)
	261	if entry is not None:
	262	for item in in_list:
	263	input_field = get_content('ns:%s' % item[0], entry)
	264	if input_field is not None:
	265	try:
	266	setattr(self, item[1], float(input_field.text.strip()))
	267	except:
	268	setattr(self, item[1], None)
	269	input_field = get_content('ns:estimate_bck', entry)
	270	if input_field is not None:
	271	try:
	272	self.estimate_bck = input_field.text.strip() == 'True'
	273	except:
	274	self.estimate_bck = False
	275
	276	# Parse inversion outputs
	277	entry = get_content('ns:outputs', node)
	278	if entry is not None:
	279	# Output parameters (scalars)
	280	for item in out_list:
	281	input_field = get_content('ns:%s' % item[0], entry)
	282	if input_field is not None:
	283	try:
	284	setattr(self, item[1], float(input_field.text.strip()))
	285	except:
	286	setattr(self, item[1], None)
	287
	288	# Look for coefficients
	289	# Format is [value, value, value, value]
	290	coeff = get_content('ns:coefficients', entry)
	291	if coeff is not None:
	292	# Remove brackets
	293	c_values = coeff.text.strip().replace('[', '')
	294	c_values = c_values.replace(']', '')
	295	toks = c_values.split()
	296	self.coefficients = []
	297	for c in toks:
	298	try:
	299	self.coefficients.append(float(c))
	300	except:
[a0e6b1b]	301	# Bad data, skip. We will count the number of
	302	# coefficients at the very end and deal with
[959eb01]	303	# inconsistencies then.
	304	pass
	305	# Sanity check
	306	if not len(self.coefficients) == self.nfunc:
	307	# Inconsistent number of coefficients.
	308	# Don't keep the data.
	309	err_msg = "InversionState.fromXML: inconsistant "
	310	err_msg += "number of coefficients: "
	311	err_msg += "%d %d" % (len(self.coefficients),
	312	self.nfunc)
	313	logger.error(err_msg)
	314	self.coefficients = None
	315
	316	# Look for covariance matrix
	317	# Format is [ [value, value], [value, value] ]
	318	coeff = get_content('ns:covariance', entry)
	319	if coeff is not None:
	320	# Parse rows
	321	rows = coeff.text.strip().split('[')
	322	self.covariance = []
	323	for row in rows:
	324	row = row.strip()
	325	if len(row) == 0: continue
	326	# Remove end bracket
	327	row = row.replace(']', '')
	328	c_values = row.split()
	329	cov_row = []
	330	for c in c_values:
	331	try:
	332	cov_row.append(float(c))
	333	except:
[a0e6b1b]	334	# Bad data, skip. We will count the number of
	335	# coefficients at the very end and deal with
[959eb01]	336	# inconsistencies then.
	337	pass
	338	# Sanity check: check the number of entries in the row
	339	if len(cov_row) == self.nfunc:
	340	self.covariance.append(cov_row)
	341	# Sanity check: check the number of rows in the covariance
	342	# matrix
	343	if not len(self.covariance) == self.nfunc:
	344	# Inconsistent dimensions of the covariance matrix.
	345	# Don't keep the data.
	346	err_msg = "InversionState.fromXML: "
	347	err_msg += "inconsistant dimensions of the "
	348	err_msg += " covariance matrix: "
	349	err_msg += "%d %d" % (len(self.covariance), self.nfunc)
	350	logger.error(err_msg)
	351	self.covariance = None
	352
	353	class Reader(CansasReader):
	354	"""
	355	Class to load a .prv P(r) inversion file
	356	"""
	357	## File type
	358	type_name = "P(r)"
	359
	360	## Wildcards
	361	type = ["P(r) files (.prv)\|.prv",
	362	"SASView files (.svs)\|.svs"]
	363	## List of allowed extensions
	364	ext = ['.prv', '.PRV', '.svs', '.SVS']
	365
	366	def __init__(self, call_back, cansas=True):
	367	"""
	368	Initialize the call-back method to be called
	369	after we load a file
	370
	371	:param call_back: call-back method
	372	:param cansas: True = files will be written/read in CanSAS format
	373	False = write CanSAS format
	374
	375	"""
	376	## Call back method to be executed after a file is read
	377	self.call_back = call_back
	378	## CanSAS format flag
	379	self.cansas = cansas
	380	self.state = None
	381
	382	def read(self, path):
	383	"""
	384	Load a new P(r) inversion state from file
	385
	386	:param path: file path
	387
	388	:return: None
	389
	390	"""
	391	if self.cansas == True:
	392	return self._read_cansas(path)
	393	else:
	394	return self._read_standalone(path)
	395
	396	def _read_standalone(self, path):
	397	"""
	398	Load a new P(r) inversion state from file.
	399	The P(r) node is assumed to be the top element.
	400
	401	:param path: file path
	402
	403	:return: None
	404
	405	"""
	406	# Read the new state from file
	407	state = InversionState()
	408	state.fromXML(file=path)
	409
	410	# Call back to post the new state
	411	self.state = state
	412	#self.call_back(state)
	413	return None
	414
	415	def _parse_prstate(self, entry):
	416	"""
	417	Read a p(r) inversion result from an XML node
	418
	419	:param entry: XML node to read from
	420
	421	:return: InversionState object
	422
	423	"""
	424	state = None
	425
	426	# Locate the P(r) node
	427	try:
	428	nodes = entry.xpath('ns:%s' % PRNODE_NAME,
	429	namespaces={'ns': CANSAS_NS})
	430	if nodes != []:
	431	# Create an empty state
	432	state = InversionState()
	433	state.fromXML(node=nodes[0])
	434	except:
	435	msg = "XML document does not contain P(r) "
	436	msg += "information.\n %s" % sys.exc_value
	437	logger.info(msg)
	438
	439	return state
	440
	441	def _read_cansas(self, path):
	442	"""
	443	Load data and P(r) information from a CanSAS XML file.
	444
	445	:param path: file path
	446
	447	:return: Data1D object if a single SASentry was found,
	448	or a list of Data1D objects if multiple entries were found,
	449	or None of nothing was found
	450
	451	:raise RuntimeError: when the file can't be opened
	452	:raise ValueError: when the length of the data vectors are inconsistent
	453
	454	"""
	455	output = []
	456
	457	if os.path.isfile(path):
	458	basename = os.path.basename(path)
	459	root, extension = os.path.splitext(basename)
	460	#TODO: eventually remove the check for .xml once
	461	# the P(r) writer/reader is truly complete.
	462	if extension.lower() in self.ext or extension.lower() == '.xml':
	463
	464	tree = etree.parse(path, parser=etree.ETCompatXMLParser())
	465	# Check the format version number
[a0e6b1b]	466	# Specifying the namespace will take care of the file
	467	#format version
[959eb01]	468	root = tree.getroot()
	469
	470	entry_list = root.xpath('/ns:SASroot/ns:SASentry',
	471	namespaces={'ns': CANSAS_NS})
	472
	473	for entry in entry_list:
	474	sas_entry, _ = self._parse_entry(entry)
	475	prstate = self._parse_prstate(entry)
	476	#prstate could be None when .svs file is loaded
	477	#in this case, skip appending to output
[7432acb]	478	if prstate is not None:
[959eb01]	479	sas_entry.meta_data['prstate'] = prstate
	480	sas_entry.filename = prstate.file
	481	output.append(sas_entry)
	482	else:
	483	raise RuntimeError, "%s is not a file" % path
	484
	485	# Return output consistent with the loader's api
	486	if len(output) == 0:
	487	return None
	488	elif len(output) == 1:
	489	# Call back to post the new state
	490	self.call_back(output[0].meta_data['prstate'], datainfo=output[0])
	491	#self.state = output[0].meta_data['prstate']
	492	return output[0]
	493	else:
	494	return output
	495
	496
	497	def write(self, filename, datainfo=None, prstate=None):
	498	"""
	499	Write the content of a Data1D as a CanSAS XML file
	500
	501	:param filename: name of the file to write
	502	:param datainfo: Data1D object
	503	:param prstate: InversionState object
	504
	505	"""
	506	# Sanity check
	507	if self.cansas == True:
	508	doc = self.write_toXML(datainfo, prstate)
	509	# Write the XML document
	510	fd = open(filename, 'w')
	511	fd.write(doc.toprettyxml())
	512	fd.close()
	513	else:
	514	prstate.toXML(file=filename)
	515
	516	def write_toXML(self, datainfo=None, state=None):
	517	"""
	518	Write toXML, a helper for write()
	519
	520	: return: xml doc
	521	"""
	522	if datainfo is None:
	523	datainfo = Data1D(x=[], y=[])
	524	elif not issubclass(datainfo.__class__, Data1D):
	525	msg = "The cansas writer expects a Data1D "
	526	msg += "instance: %s" % str(datainfo.__class__.__name__)
	527	raise RuntimeError, msg
	528
	529	# Create basic XML document
	530	doc, sasentry = self._to_xml_doc(datainfo)
	531
	532	# Add the invariant information to the XML document
	533	if state is not None:
	534	doc = state.toXML(doc=doc, entry_node=sasentry)
	535
	536	return doc

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

SasView

source: sasview/src/sas/sasgui/perspectives/pr/inversion_state.py @ 3a08f369

Download in other formats: