maped_average.py 2.58 KB
title = "MapAverage"
tip = "average maps with respect to the same center"
onein = False

import numpy as np
from nodmath import extract
from scipy import stats

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 _

def nint(x):
    if x > 0: return int(x+0.5)
    else: return int(x-0.5)


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):
            posx = IntItem('center x')
            posy = IntItem('center y')
            naxis1 = IntItem('shape x')
            naxis2 = IntItem('shape y')
        name = title.replace(" ", "")
        if args == {}:
           param = FuncParam(_(title), "Select center pixel and dimensions:\n")
        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 function(self, ms, p):
        px = []
        py = []
        nx = []
        ny = []
        i = 0
        j = 0
        cpix1 = None
        cpix2 = None
        for m in ms:
            m = self.parent.removeNANedges(m)
            if cpix1 == None:
               cpix1 = m.header['CRPIX1']
            if cpix2 == None:
               cpix2 = m.header['CRPIX2']
            cpix1 = max(cpix1, m.header['CRPIX1'])
            cpix2 = max(cpix2, m.header['CRPIX2'])
            pix = nint(2*m.header['CRPIX1'])/2.0
            piy = nint(2*m.header['CRPIX2'])/2.0
            if pix*abs(m.header['CDELT1']) > 180.0:
               pix -= abs(360.0/m.header['CDELT1'])
            px.append(pix)
            py.append(piy)
            nx.append(2*max(abs(pix - m.header['NAXIS1']), abs(pix)))
            ny.append(2*max(abs(m.header['NAXIS2'] - piy), abs(piy)))   
        cols = int(max(nx))
        rows = int(max(ny))
        maps = []
        m = ms[0]
        for n in range(len(px)):
            d = extract(ms[n].data, shape=(rows, cols), position=(py[n], px[n]))
            #ms[n].data = d
            maps.append(d)
        m.data = stats.nanmean(maps, axis=0)
        m = self.parent.removeNANedges(m)
        m.header['CRPIX1'] = cpix1
        m.header['CRPIX2'] = cpix2
        return m, p