# -*- coding: utf-8 -*-
"""
/***************************************************************************
 BeachBall
                                 A QGIS plugin
 This plugin creates beachball to represent focal mechanisms from strike, dip, and rake.
 Generated by Plugin Builder: http://g-sherman.github.io/Qgis-Plugin-Builder/
                              -------------------
        begin                : 2025-04-01
        git sha              : $Format:%H$
        copyright            : (C) 2025 by GeoSeisUtilities
        email                : zazo.dt@gmail.com
 ***************************************************************************/

/***************************************************************************
 *                                                                         *
 *   This program is free software; you can redistribute it and/or modify  *
 *   it under the terms of the GNU General Public License as published by  *
 *   the Free Software Foundation; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 ***************************************************************************/
"""
from qgis.core import QgsProject, QgsField, Qgis, QgsSvgMarkerSymbolLayer, QgsMarkerSymbol, QgsSymbolLayer, QgsProperty
from qgis.PyQt.QtCore import QSettings, QTranslator, QCoreApplication, QVariant
from qgis.PyQt.QtGui import QIcon, QColor
from qgis.PyQt.QtWidgets import QAction, QFileDialog, QMessageBox
from qgis.core import QgsVectorLayer
import processing

# Initialize Qt resources from file resources.py
from .resources import *
# Import the code for the dialog
from .beach_ball_dialog import BeachBallDialog
import os.path
import sys


class BeachBall:
    """QGIS Plugin Implementation."""

    def __init__(self, iface):
        """Constructor.

        :param iface: An interface instance that will be passed to this class
            which provides the hook by which you can manipulate the QGIS
            application at run time.
        :type iface: QgsInterface
        """
        # Save reference to the QGIS interface
        self.iface = iface
        # initialize plugin directory
        self.plugin_dir = os.path.dirname(__file__)
        # initialize libs directory
        self.libs_path = os.path.join(self.plugin_dir, 'libs')
        if self.libs_path not in sys.path:
            sys.path.insert(0, self.libs_path)
        # initialize locale
        locale = QSettings().value('locale/userLocale')[0:2]
        locale_path = os.path.join(
            self.plugin_dir,
            'i18n',
            'BeachBall_{}.qm'.format(locale))

        if os.path.exists(locale_path):
            self.translator = QTranslator()
            self.translator.load(locale_path)
            QCoreApplication.installTranslator(self.translator)

        # Declare instance attributes
        self.actions = []
        self.menu = self.tr(u'&BeachBall')

        # Check if plugin was started the first time in current QGIS session
        # Must be set in initGui() to survive plugin reloads
        self.first_start = None

    # noinspection PyMethodMayBeStatic
    def tr(self, message):
        """Get the translation for a string using Qt translation API.

        We implement this ourselves since we do not inherit QObject.

        :param message: String for translation.
        :type message: str, QString

        :returns: Translated version of message.
        :rtype: QString
        """
        # noinspection PyTypeChecker,PyArgumentList,PyCallByClass
        return QCoreApplication.translate('BeachBall', message)


    def add_action(
        self,
        icon_path,
        text,
        callback,
        enabled_flag=True,
        add_to_menu=True,
        add_to_toolbar=True,
        status_tip=None,
        whats_this=None,
        parent=None):
        """Add a toolbar icon to the toolbar.

        :param icon_path: Path to the icon for this action. Can be a resource
            path (e.g. ':/plugins/foo/bar.png') or a normal file system path.
        :type icon_path: str

        :param text: Text that should be shown in menu items for this action.
        :type text: str

        :param callback: Function to be called when the action is triggered.
        :type callback: function

        :param enabled_flag: A flag indicating if the action should be enabled
            by default. Defaults to True.
        :type enabled_flag: bool

        :param add_to_menu: Flag indicating whether the action should also
            be added to the menu. Defaults to True.
        :type add_to_menu: bool

        :param add_to_toolbar: Flag indicating whether the action should also
            be added to the toolbar. Defaults to True.
        :type add_to_toolbar: bool

        :param status_tip: Optional text to show in a popup when mouse pointer
            hovers over the action.
        :type status_tip: str

        :param parent: Parent widget for the new action. Defaults None.
        :type parent: QWidget

        :param whats_this: Optional text to show in the status bar when the
            mouse pointer hovers over the action.

        :returns: The action that was created. Note that the action is also
            added to self.actions list.
        :rtype: QAction
        """

        icon = QIcon(icon_path)
        action = QAction(icon, text, parent)
        action.triggered.connect(callback)
        action.setEnabled(enabled_flag)

        if status_tip is not None:
            action.setStatusTip(status_tip)

        if whats_this is not None:
            action.setWhatsThis(whats_this)

        if add_to_toolbar:
            # Adds plugin icon to Plugins toolbar
            self.iface.addToolBarIcon(action)

        if add_to_menu:
            self.iface.addPluginToMenu(
                self.menu,
                action)

        self.actions.append(action)

        return action

    def initGui(self):
        """Create the menu entries and toolbar icons inside the QGIS GUI."""

        icon_path = ':/plugins/beach_ball/icon.png'
        self.add_action(
            icon_path,
            text=self.tr(u'Draw beachball for focal mechanisms'),
            callback=self.run,
            parent=self.iface.mainWindow())

        # will be set False in run()
        self.first_start = True

    def unload(self):
        """Removes the plugin menu item and icon from QGIS GUI."""
        for action in self.actions:
            self.iface.removePluginMenu(
                self.tr(u'&BeachBall'),
                action)
            self.iface.removeToolBarIcon(action)

    def select_output_file(self):
        folder = QFileDialog.getExistingDirectory(
            self.dlg)
        self.dlg.lineEdit.setText(folder)

    def getOutput(self):
        # Return the output file path, or use a temporary file if blank
        output = self.dlg.lineEdit.text().strip()
        if output == '':
            output = '/'.join(QgsProject.instance().fileName().split('/')[:-1])
            if output == '':
                output='/home'
        output = output.replace('\\', '/')
        return output

    def run(self):
        """Run method that performs all the real work"""
        
        # Create the dialog with elements (after translation) and keep reference
        # Only create GUI ONCE in callback, so that it will only load when the plugin is started
        if self.first_start == True:
            self.first_start = False
            self.dlg = BeachBallDialog()
            self.dlg.pushButton.clicked.connect(self.select_output_file)
            # Set color of the color buttons
            self.dlg.colorButton_2.setColor(QColor(255, 0, 0))
            self.dlg.colorButton_3.setColor(QColor(255, 153, 51))
            self.dlg.colorButton_4.setColor(QColor(0, 204, 0))
            self.dlg.colorButton_5.setColor(QColor(0, 255, 255))
            self.dlg.colorButton_6.setColor(QColor(0, 0, 255))

        def showError(message):
            QMessageBox.critical(None, "Error", message)

        def get_selected_ly():
            # Clear the contents of the comboBox from previous runs
            self.dlg.comboBox_2.clear()
            self.dlg.comboBox_3.clear()
            self.dlg.comboBox_4.clear()
            self.dlg.comboBox_5.clear()
            self.dlg.comboBox_6.clear()
            # Get layer fields
            layer = self.dlg.comboBox.currentData()
            if layer == None: # Invalid layer
                return
            if not isinstance(layer, QgsVectorLayer):
                return
            fields = layer.fields().names()
            # Populate the comboBox with names of the fields
            self.dlg.comboBox_2.addItems(fields)
            self.dlg.comboBox_3.addItems(fields)
            self.dlg.comboBox_4.addItems(fields)
            self.dlg.comboBox_5.addItems(fields)
            self.dlg.comboBox_6.addItem("")
            self.dlg.comboBox_6.addItems(fields)
            self.dlg.comboBox_6.setCurrentIndex(0)

        # Needed modules
        import os
        import shutil
        from qgis.PyQt.QtWidgets import QDialogButtonBox
        import numpy as np
        from obspy.imaging.beachball import beachball

        # Fetch the currently loaded layers
        layers = QgsProject.instance().mapLayers().values()
        # Clear the contents of the comboBox from previous runs
        self.dlg.comboBox.clear()
        # Populate the comboBox with names of all the loaded layers
        for layer in layers:
            self.dlg.comboBox.addItem(layer.name(), layer)
        get_selected_ly()
        self.dlg.comboBox.currentIndexChanged.connect(get_selected_ly)

        def execute_processing():
            output_path = self.getOutput()
            if os.path.isdir(output_path) == False:
                showError(f"Output path is invalid:\n Path:{output_path}")

            layer = self.dlg.comboBox.currentData()
            id_idx = self.dlg.comboBox_2.currentText()
            strike_idx = self.dlg.comboBox_3.currentText()
            dip_idx = self.dlg.comboBox_4.currentText()
            rake_idx = self.dlg.comboBox_5.currentText()
            mag_idx = self.dlg.comboBox_6.currentText()
            layer.selectAll()
            clone_layer = processing.run("native:saveselectedfeatures", {'INPUT': layer, 'OUTPUT': 'memory:'})['OUTPUT']
            layer.removeSelection()
            QgsProject.instance().addMapLayer(clone_layer)
            NF_col = self.dlg.colorButton_2.color().name()
            NS_col = self.dlg.colorButton_3.color().name()
            SS_col = self.dlg.colorButton_4.color().name()
            TS_col = self.dlg.colorButton_5.color().name()
            TF_col = self.dlg.colorButton_6.color().name()
            if os.path.exists(output_path+'/bb_svg'):
                shutil.rmtree(output_path+'/bb_svg' + '/')
            os.mkdir(output_path+'/bb_svg')
            # Start editing
            clone_layer.startEditing()
            field_name = "Kin"
            while clone_layer.dataProvider().fieldNameIndex(field_name) != -1:
                field_name = field_name +'_'
            clone_layer.dataProvider().addAttributes([QgsField(field_name, QVariant.String)])
            clone_layer.updateFields()
            id_new_col= clone_layer.dataProvider().fieldNameIndex(field_name)
            field_name = "Ptrend"
            while clone_layer.dataProvider().fieldNameIndex(field_name) != -1:
                field_name = field_name +'_'
            clone_layer.dataProvider().addAttributes([QgsField(field_name, QVariant.Double)])
            clone_layer.updateFields()
            id_Ptrend= clone_layer.dataProvider().fieldNameIndex(field_name)
            field_name = "Pplunge"
            while clone_layer.dataProvider().fieldNameIndex(field_name) != -1:
                field_name = field_name +'_'
            clone_layer.dataProvider().addAttributes([QgsField(field_name, QVariant.Double)])
            clone_layer.updateFields()
            id_Pplunge= clone_layer.dataProvider().fieldNameIndex(field_name)
            field_name = "Ttrend"
            while clone_layer.dataProvider().fieldNameIndex(field_name) != -1:
                field_name = field_name +'_'
            clone_layer.dataProvider().addAttributes([QgsField(field_name, QVariant.Double)])
            clone_layer.updateFields()
            id_Ttrend= clone_layer.dataProvider().fieldNameIndex(field_name)
            field_name = "Tplunge"
            while clone_layer.dataProvider().fieldNameIndex(field_name) != -1:
                field_name = field_name +'_'
            clone_layer.dataProvider().addAttributes([QgsField(field_name, QVariant.Double)])
            clone_layer.updateFields()
            id_Tplunge= clone_layer.dataProvider().fieldNameIndex(field_name)
            field_name = "Btrend"
            while clone_layer.dataProvider().fieldNameIndex(field_name) != -1:
                field_name = field_name +'_'
            clone_layer.dataProvider().addAttributes([QgsField(field_name, QVariant.Double)])
            clone_layer.updateFields()
            id_Btrend= clone_layer.dataProvider().fieldNameIndex(field_name)
            field_name = "Bplunge"
            while clone_layer.dataProvider().fieldNameIndex(field_name) != -1:
                field_name = field_name +'_'
            clone_layer.dataProvider().addAttributes([QgsField(field_name, QVariant.Double)])
            clone_layer.updateFields()
            id_Bplunge= clone_layer.dataProvider().fieldNameIndex(field_name)
            field_name = "SHmax"
            while clone_layer.dataProvider().fieldNameIndex(field_name) != -1:
                field_name = field_name +'_'
            clone_layer.dataProvider().addAttributes([QgsField(field_name, QVariant.Double)])
            clone_layer.updateFields()
            id_SHmax= clone_layer.dataProvider().fieldNameIndex(field_name)
            field_name = "SVGpath"
            while clone_layer.dataProvider().fieldNameIndex(field_name) != -1:
                field_name = field_name +'_'
            clone_layer.dataProvider().addAttributes([QgsField(field_name, QVariant.String)])
            clone_layer.updateFields()
            id_svg_col= clone_layer.dataProvider().fieldNameIndex(field_name)
            # Calculating beachball
            nvi = 0
            for feature in clone_layer.getFeatures():
                try:
                    id_v = feature[id_idx]
                except:
                    svg_file = f'{output_path}/bb_svg/NotValid{nvi}.svg'
                    clone_layer.changeAttributeValue(feature.id(), id_svg_col, svg_file)
                    beachball([0, 0, 0], size=0, linewidth=0, facecolor='white', outfile=svg_file)
                    nvi += 0
                    continue
                try:
                    strike = float(feature[strike_idx])
                except:
                    svg_file = f'{output_path}/bb_svg/NotValid{nvi}.svg'
                    clone_layer.changeAttributeValue(feature.id(), id_svg_col, svg_file)
                    beachball([0, 0, 0], size=0, linewidth=0, facecolor='white', outfile=svg_file)
                    nvi += 0
                    continue
                try:
                    dip = float(feature[dip_idx])
                except:
                    svg_file = f'{output_path}/bb_svg/NotValid{nvi}.svg'
                    clone_layer.changeAttributeValue(feature.id(), id_svg_col, svg_file)
                    beachball([0, 0, 0], size=0, linewidth=0, facecolor='white', outfile=svg_file)
                    nvi += 0
                    continue
                try:
                    rake = float(feature[rake_idx])
                except:
                    svg_file = f'{output_path}/bb_svg/NotValid{nvi}.svg'
                    clone_layer.changeAttributeValue(feature.id(), id_svg_col, svg_file)
                    beachball([0, 0, 0], size=0, linewidth=0, facecolor='white', outfile=svg_file)
                    nvi += 0
                    continue
                n1 = -np.sin(np.radians(dip)) * np.sin(np.radians(strike))
                n2 = np.sin(np.radians(dip)) * np.cos(np.radians(strike))
                n3 = -np.cos(np.radians(dip))
                n = np.array([n1, n2, n3])
                u1 = (np.cos(np.radians(rake)) * np.cos(np.radians(strike)) +
                    np.cos(np.radians(dip)) * np.sin(np.radians(rake)) * np.sin(np.radians(strike)))
                u2 = (np.cos(np.radians(rake)) * np.sin(np.radians(strike)) -
                    np.cos(np.radians(dip)) * np.sin(np.radians(rake)) * np.cos(np.radians(strike)))
                u3 = -np.sin(np.radians(rake)) * np.sin(np.radians(dip))
                u = np.array([u1, u2, u3])
                P_osa = np.subtract(n, u) / np.linalg.norm(np.subtract(n, u))
                T_osa = np.add(n, u) / np.linalg.norm(np.add(n, u))
                B_osa = np.cross(T_osa, P_osa) / np.linalg.norm(np.cross(T_osa, P_osa))
                if P_osa[2] > 0:
                    P_osa = -P_osa
                if T_osa[2] > 0:
                    T_osa = -T_osa
                if B_osa[2] > 0:
                    B_osa = -B_osa
                if P_osa[2] > 0:
                    correctP = 0
                else:
                    correctP = 180
                if T_osa[2] > 0:
                    correctT = 0
                else:
                    correctT = 180
                if B_osa[2] > 0:
                    correctB = 0
                else:
                    correctB = 180
                P_theta = np.degrees(np.arcsin(abs(P_osa[2])))
                T_theta = np.degrees(np.arcsin(abs(T_osa[2])))
                B_theta = np.degrees(np.arcsin(abs(B_osa[2])))
                P_azimuth = (correctP + (np.degrees(np.arctan2((P_osa[1]), P_osa[0])) + 360.0)) % 360.0
                T_azimuth = (correctT + (np.degrees(np.arctan2(T_osa[1], T_osa[0])) + 360.0)) % 360.0
                B_azimuth = (correctB + (np.degrees(np.arctan2((B_osa[1]), B_osa[0])) + 360.0)) % 360.0
                SH_p = P_osa
                SH_p[2] = 0
                SH_t = T_osa
                SH_t[2] = 0
                if np.linalg.norm(SH_p) >= np.linalg.norm(SH_t):
                    SH = SH_p
                else:
                    SH = SH_t
                SH /= np.linalg.norm(SH)
                y, x, _ = SH / np.linalg.norm(SH)
                SHmax = np.degrees(np.arctan2(x, y)) % 360
                if P_theta>=52 and T_theta<=35:
                    kin='NF'
                    col = NF_col
                elif P_theta>=40 and P_theta<52 and T_theta<=20:
                    kin='NS'
                    col = NS_col
                elif P_theta<40 and T_theta<=20 and B_theta>=45:
                    kin='SS'
                    col = SS_col
                elif P_theta<=20 and T_theta<40 and B_theta>=45:
                    kin='SS'
                    col = SS_col
                elif P_theta<=20 and T_theta>=40 and T_theta<52:
                    kin='TS'
                    col = TS_col
                elif P_theta<=35 and T_theta>=52:
                    kin='TF'
                    col = TF_col
                else:
                    kin='UK'
                    col = 'grey'
                # Update new field
                for val in [np.nan, "", None, np.inf]:
                    if P_azimuth == val or float(P_azimuth) == val:
                        P_azimuth = None
                    else:
                        P_azimuth = float(P_azimuth)
                    if P_theta == val or float(P_theta) == val:
                        P_theta = None
                    else:
                        P_theta = float(P_theta)
                    if T_azimuth == val or float(T_azimuth) == val:
                        T_azimuth = None
                    else:
                        T_azimuth = float(T_azimuth)
                    if T_theta == val or float(T_theta) == val:
                        T_theta = None
                    else:
                        T_theta = float(T_theta)
                    if B_azimuth == val or float(B_azimuth) == val:
                        B_azimuth = None
                    else:
                        B_azimuth = float(B_azimuth)
                    if B_theta == val or float(B_theta) == val:
                        B_theta = None
                    else:
                        B_theta = float(B_theta)
                    if SHmax == val or float(SHmax) == val:
                        SHmax = None
                    else:
                        SHmax = float(SHmax)
                clone_layer.changeAttributeValue(feature.id(), id_new_col, kin)
                clone_layer.changeAttributeValue(feature.id(), id_Ptrend, round(P_azimuth,2))
                clone_layer.changeAttributeValue(feature.id(), id_Pplunge, round(P_theta,2))
                clone_layer.changeAttributeValue(feature.id(), id_Ttrend, round(T_azimuth,2))
                clone_layer.changeAttributeValue(feature.id(), id_Tplunge, round(T_theta,2))
                clone_layer.changeAttributeValue(feature.id(), id_Btrend, round(B_azimuth,2))
                clone_layer.changeAttributeValue(feature.id(), id_Bplunge, round(B_theta,2))
                clone_layer.changeAttributeValue(feature.id(), id_SHmax, round(SHmax,2))
                # Plotting beachball
                svg_file = f'{output_path}/bb_svg/{id_v}.svg'
                clone_layer.changeAttributeValue(feature.id(), id_svg_col, svg_file)
                beachball([strike, dip, rake], size=200, linewidth=2, facecolor=col, outfile=svg_file)

            clone_layer.commitChanges()
            # Adding custom style
            symbol = QgsMarkerSymbol.createSimple({})
            svgLayer = QgsSvgMarkerSymbolLayer(svg_file)
            svgLayer.setDataDefinedProperty(QgsSymbolLayer.PropertyName, QgsProperty.fromExpression('"SVGpath"'))
            if mag_idx and mag_idx.strip():
                exp_string = 'coalesce("{0}" * 2, 7)'.format(mag_idx)
                svgLayer.setDataDefinedProperty(QgsSymbolLayer.PropertySize,
                    QgsProperty.fromExpression(exp_string))
            else:
                svgLayer.setSize(7.0)
            symbol.changeSymbolLayer(0, svgLayer)
            clone_layer.renderer().setSymbol(symbol)
            clone_layer.triggerRepaint()
            self.iface.messageBar().pushMessage("BeachBall computed",
                level=Qgis.Success, duration=3)
            self.dlg.accept()

        ok_btn = self.dlg.button_box.button(QDialogButtonBox.Ok)
        try:
            ok_btn.clicked.disconnect()
        except:
            pass
        ok_btn.clicked.connect(execute_processing)
        self.dlg.exec_()
