source: sasview/prview/perspectives/pr/inversion_state.py @ 6ab0ad1

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 6ab0ad1 was b1f7ec6, checked in by Mathieu Doucet <doucetm@…>, 15 years ago

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