source: sasview/src/sas/sasgui/perspectives/pr/inversion_state.py @ 25dd9c9

Last change on this file since 25dd9c9 was 7432acb, checked in by andyfaff, 8 years ago

MAINT: search+replace '!= None' by 'is not None'

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