"""
Geosoft grid and image utilities.

.. seealso:: `geosoft.gxpy.grid`, `geosoft.gxapi.GXIMG`, `geosoft.gxapi.GXIMU`

.. note::

    Regression tests provide usage examples:     
    `Tests <https://github.com/GeosoftInc/gxpy/blob/master/geosoft/gxpy/tests/test_grid_utility.py>`_
    
"""
import os
import numpy as np
import math

import geosoft
import geosoft.gxapi as gxapi
from . import gx as gx
from . import vv as gxvv
from . import grid as gxgrd
from . import map as gxmap
from . import view as gxv
from . import group as gxgrp
from . import gdb as gxgdb
from . import geometry as gxgeo
from . import utility as gxu
from . import geometry_utility as gxgeou
from . import grid_fft as gxfft

__version__ = geosoft.__version__

BOOL_AND = gxapi.IMU_BOOL_OPT_AND
BOOL_OR = gxapi.IMU_BOOL_OPT_OR
BOOL_XOR = gxapi.IMU_BOOL_OPT_XOR
BOOL_SIZE_GRID1 = gxapi.IMU_BOOL_SIZING_0
BOOL_SIZE_GRID2 = gxapi.IMU_BOOL_SIZING_1
BOOL_SIZE_MIN = gxapi.IMU_BOOL_SIZING_MIN
BOOL_SIZE_MAX = gxapi.IMU_BOOL_SIZING_MAX
BOOL_OVERLAP_AVERAGE = gxapi.IMU_BOOL_OLAP_AVE
BOOL_OVERLAP_GRID1 = gxapi.IMU_BOOL_OLAP_1
BOOL_OVERLAP_GRID2 = gxapi.IMU_BOOL_OLAP_2

TREND_EDGE = gxapi.IMU_TREND_EDGE
TREND_ALL = gxapi.IMU_TREND_ALL
DERIVATIVE_X = 0
DERIVATIVE_Y = 1
DERIVATIVE_Z = 2
DERIVATIVE_XY = 3
DERIVATIVE_XYZ = 4
TILT_ANGLE = 5
RETURN_PPOINT = 0
RETURN_LIST_OF_PPOINT = 1
RETURN_GDB = 2


def _t(s):
    return geosoft.gxpy.system.translate(s)


class GridUtilityException(geosoft.GXRuntimeError):
    """
    Exceptions from `geosoft.gxpy.grid_utility`.

    .. versionadded:: 9.4
    """
    pass


def remove_trend(grid, file_name=None, method=TREND_EDGE, overwrite=False):
    """
    Calculate a polynomial trend surface and return trend-removed grid.

    :param grid:        `geosoft.gxpy.grid.Grid` instance, or a file name
    :param file_name:   trend-removed grid file name, if `None` a temporary grid is created.
    :param method:      base trend on `TREND_EDGE` for edge data or `TREND_ALL` for all data
    :param overwrite:   True to overwrite existing file_name
    :return:            `Grid` instance
    
    .. versionadded 9.4
    """

    if not isinstance(grid, gxgrd.Grid):
        grid = gxgrd.Grid.open(grid, dtype=np.float64, mode=gxgrd.FILE_READ)

    # need GS_DOUBLE grids
    if grid.gxtype != gxapi.GS_DOUBLE:
        ing = grid.copy(grid, gx.gx().temp_file('.grd(GRD)'), dtype=np.float64)
        ing.delete_files()
    else:
        ing = grid

    if file_name is None:
        file_name = gx.gx().temp_file('.grd(GRD)')
    dtg = grid.new(file_name=file_name, properties=ing.properties(), overwrite=overwrite)
    gxapi.GXIMU.grid_trnd(ing.gximg, dtg.gximg, 0, method, 1,
                          gxapi.GXVM.create(gxapi.GS_REAL, 10), 3)
    return dtg


def derivative(grid, derivative_type, file_name=None, overwrite=False, dtype=None, fft=True):
    """
    Return a derivative of a grid or a tilt-angle

    :param grid:            `geosoft.gxpy.grid.Grid` instance, or a file name
    :param derivative_type:  Which derivative to calculate, (grid_data_uom/distance):
    
        ========================== ====================================================================
        DERIVATIVE_X               in the grid X direction
        DERIVATIVE_Y               in the grid Y direction
        DERIVATIVE_Z               in the grid Z direction
        DERIVATIVE_XYZ             the total derivative sqrt(dx**2 + dy**2 + dz**2) (analytic signal)
        TILT_ANGLE                 tilt angle, atan2(dz, sqrt(dx**2 + dy**2)) (radians)
        ========================== ====================================================================
    
    :param file_name:   returned derivative file name, `None` for a temporary file
    :param overwrite:   True to overwrite existing file
    :param dtype:       dtype for the return grid, default is the same as the passed grid.
    :param fft:         `False` calculate Z derivative with a space-domain convolution rather than an FFT.
    :return:            `geosoft.gxpy.grid.Grid` instance that contains the derivative result

    .. note:: Derivative units_of_measure are grid_unit_of_measure / distance, except for the tilt angle,
        which is radians.

        Horizontal derivatives are calculated in the space domain based on the difference between
        neighboring cell values, and the Z derivative can be calculated using an FFT or a 5x5 space-domain
        convolution.  An FFT calculation will generally produce a better result and it will be able to
        work with longer wavelengths, but at the expense of speed and edge effects in cases of very powerful
        anomalies along the edge of a grid.

    .. versionadded 9.4
    """

    def vertical_derivative(g, dt):

        # float64 grids for grid_vd
        if not isinstance(g, gxgrd.Grid):
            g = gxgrd.Grid.open(g, dtype=np.float64, mode=gxgrd.FILE_READ)

        if g.dtype != np.float64:
            g = g.copy(g, gx.gx().temp_file('.grd(GRD)'), dtype=np.float32, overwrite=True)
            g.delete_files()

        if fft:
            with gxfft.GridFFT(g) as gfft:
                gfft.filter(filters=['DRVZ 1'])