misc_roi.py 3.76 KB
title = "Statistics"
tip = "calculates statistical value from map extraction"
onein = True

import numpy as np

from guidata.qt.QtGui import QMessageBox
from guidata.dataset.datatypes import DataSet
from guidata.dataset.dataitems import (IntItem, StringItem, ChoiceItem, FloatItem, BoolItem)
from guiqwt.config import _

class NOD3_App:

    def __init__(self, parent):
        self.parent = parent
        self.parent.activateWindow()

    def Error(self, msg):
        QMessageBox.critical(self.parent.parent(), title,
                              _(u"Error:")+"\n%s" % str(msg))

    def compute_app(self, **args):
        class FuncParam(DataSet):
            row1 = IntItem(_("First row index"), default=0, min=-1)
            row2 = IntItem(_("Last row index"), default=-1, min=-1)
            col1 = IntItem(_("First column index"), default=0, min=-1)
            col2 = IntItem(_("Last column index"), default=-1, min=-1)
        name = title.replace(" ", "")
        if args == {}:
           param = FuncParam(_(title), "select ROI from map")
        else:
           param = self.parent.ScriptParameter(name, args)

        # if no parameter needed set param to None. activate next line
        #param = None
        self.parent.compute_11(name, lambda m, p: self.function(m, p), param, onein) 

    def PrintOut(self, data, filename='Statistics.out'):
        try:
           f = open(filename, 'w')
        except:
           self.Error("no permisson to write on disk")
           return
        f.write('# Statistic output:\n')
        f.write('\n')
        f.write(str('xmin = %-6d      xmax = %-6d\n' % (int(data[0]), int(data[1]))))
        f.write(str('ymin = %-6d      ymax = %-6d\n' % (int(data[2]), int(data[3]))))
        f.write(str('min = %f  max = %f\n' % (float(data[4]), float(data[5]))))
        f.write(str('mean = %f\n' % float(data[6])))
        f.write(str('median = %f\n' % float(data[7])))
        f.write(str('sum = %f\n' % float(data[8])))
        if type(data[9]) == str:
           f.write(str('integral = %s\n' % data[9]))
        else:
           f.write(str('integral = %f\n' % float(data[9])))
        f.write(str('variance = %f\n' % float(data[10])))
        f.write(str('rms = %f\n' % float(data[11])))
        f.close()

    def function(self, m, p):
        rows, cols = m.data.shape
        if p.row2 == -1: p.row2 = rows
        if p.col2 == -1: p.col2 = cols
        roi = m.data[p.row1:p.row2, p.col1:p.col2]
        dmin = np.nanmin(roi)
        dmax = np.nanmax(roi)
        mask = ~np.isnan(roi)
        dmean = np.mean(roi[mask])
        dmed = np.median(roi[mask])
        dvar = np.var(roi[mask], ddof=1)
        drms = np.std(roi[mask], ddof=1)
        dsum = np.sum(roi[mask])
        try:
           #beam = 2*np.pi/(8*np.log(2.0)) * (m.header['BMAJ'] * m.header['BMIN']) / abs(m.header['CDELT1'] * m.header['CDELT2']) 
           beam = abs(m.header['CDELT1'] * m.header['CDELT2']) / (m.header['BMAJ'] * m.header['BMIN']) / 1.22**2
           #integral = dsum/beam
           integral = dsum*beam
           integral = str("%.5f" % integral)
        except:
           integral = "no beam-size available"
        QMessageBox.about(self.parent.parent(), title,
            """<b>xmin = </b>%d <b>xmax = </b>%d
               <br><b>ymin = </b>%d <b>ymax = </b>%d
               <br><b>min = </b>%.5f <br><b>max = </b>%.5f
               <br><b>mean = </b>%.5f <br><b>median = </b>%.5f
               <br><b>sum = </b>%.5f <br><b>integral = </b>%s
               <br><b>variance = </b>%.5f
               <br><b>rms = </b>%.5f""" % \
               (p.col1, p.col2, p.row1, p.row2, dmin, dmax, dmean, dmed, dsum, integral, dvar, drms))
        inp = (p.col1, p.col2, p.row1, p.row2, dmin, dmax, dmean, dmed, dsum, integral, dvar, drms)
        self.PrintOut(inp)
        return [], None