Extractor_dialog.py

# -*- coding: utf-8 -*-
"""
/***************************************************************************
 ExtractorDialog
                                 A QGIS plugin
 Extraction tool to extract values from raster in many way
 Generated by Plugin Builder: http://g-sherman.github.io/Qgis-Plugin-Builder/
                             -------------------
        begin                : 2023-03-08
        git sha              : $Format:%H$
        copyright            : (C) 2023 by Stanislas Mahussi GANDAHO
        email                : gandahostanmah@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.                                   *
 *                                                                         *
 ***************************************************************************/
"""

import os
from qgis.PyQt import uic, QtWidgets
## my own import
from qgis.gui import QgsFileWidget, QgsCheckableComboBox
from qgis.core import QgsVectorLayer, QgsRasterLayer, QgsProject
from qgis.utils import iface

from PyQt5.QtCore import Qt
from osgeo import gdal, ogr, osr

import pandas as pd
import numpy as np
from .src.Zonal_stats import zonal_data, groupby_agg
from .src.AddDelim import Delim
from .src.DialogMessage import (missed_raster_vector, missed_raster_path, missed_vector_path, Unoverlable,
not_point, not_polygon, not_line, empty_folder, no_folder, is_directory, no_stat)
# Declare file path and extension
file_path = folder_path = ""
rast_ext = ('.tif', '.TIF', '.TIFF', '.tiff')
### Pandas model ###
from .src.PandasModel import PandasModel

# Loads .ui file that PyQt can populate from Qt Designer
FORM_CLASS, _ = uic.loadUiType(os.path.join(
    os.path.dirname(__file__), 'Extractor_dialog_base.ui'))

class ExtractorDialog(QtWidgets.QDialog, FORM_CLASS):
    def __init__(self, parent=None):
        """Constructor."""
        super(ExtractorDialog, self).__init__(parent)
        # Set up the user interface from Designer through FORM_CLASS.
        # After self.setupUi() you can access any designer object by doing
        # self.<objectname>, and you can use autoconnect slots - see
        # http://qt-project.org/doc/qt-4.8/designer-using-a-ui-file.html
        # #widgets-and-dialogs-with-auto-connect
        self.setupUi(self)
        self.setWindowFlags(Qt.WindowMinimizeButtonHint | Qt.WindowMaximizeButtonHint | Qt.WindowCloseButtonHint)
        qgis_font = iface.mainWindow().font() 
        self.setFont(qgis_font)
        ## set diferent path acess
        self.SingleRasterCheck.stateChanged.connect(lambda: self.PathToFile.setEnabled(self.SingleRasterCheck.isChecked()))
        self.MultipleRasterCheck.stateChanged.connect(lambda: self.PathToFile.setEnabled(self.MultipleRasterCheck.isChecked()))
        self.SingleRasterCheck.stateChanged.connect(self.path_to_file)
        self.MultipleRasterCheck.stateChanged.connect(self.path_to_file)
        self.PathToFile = self.findChild(QgsFileWidget, 'PathToFile')
        self.PathToFile.fileChanged.connect(lambda: self.path_to_file(self.PathToFile))

        # Create the table view based on QTableView widget named RasterTableView
        self.raster_table_view = self.findChild(QtWidgets.QTableView, 'RasterTableView')
        # Create the extract button
        self.Extract.clicked.connect(self.display_raster_data)
        self.Extract.clicked.connect(self.path_to_file)
        self.Extract.clicked.connect(self.handle_shp_feature_type_error)

        ## CLEAR TABLE VIEW
        self.Clear = self.findChild(QtWidgets.QPushButton, 'ClearRasterTable')
        self.Clear.clicked.connect(self.clear_raster_table)
        # Store the original data in a variable
        self.originalDataFrame = pd.DataFrame()
        # Load the original data into the table view when the application starts
        model = PandasModel(self.originalDataFrame)
        self.raster_table_view.setModel(model)
        ## SAVE TABLE
        self.Save = self.findChild(QtWidgets.QPushButton, 'SaveRasterTable')
        self.Save.clicked.connect(self.save_raster_table)
        self.Extract.clicked.connect(self.enableSave)
        self.Clear.clicked.connect(self.unableSave)

        ## VECTOR FILE CHOICE
        self.PathToVectorFile = self.findChild(QgsFileWidget, 'PathToVectorFile')
        self.PathToVectorFile.setVisible(False)
        self.VectorType.currentIndexChanged.connect(self.selectVectorFile)
        #
        self.PathToVectorFile.setStorageMode(QgsFileWidget.GetFile)
        self.PathToVectorFile.fileChanged.connect(lambda: self.vector_file_path(self.PathToVectorFile))
        self.PathToVectorFile.fileChanged.connect(self.field_name)
        ## VECTOR FIELD CHOICE
        self.FieldName = self.findChild(QtWidgets.QComboBox, 'FieldName')
        self.FieldName.setVisible(False)
        self.FieldName.currentIndexChanged.connect(self.field_chosen)
        ## WHOLE DATA OR STAT
        self.WholeOrStat = self.findChild(QtWidgets.QComboBox, 'WholeOrStat')
        self.WholeOrStat.setVisible(False)
        self.VectorType.currentIndexChanged.connect(self.whole_or_stat)
        self.MultipleRasterCheck.stateChanged.connect(self.remove_whole_or_stat)
        ## STATS CHOSEN
        StatsChoice = self.findChild(QgsCheckableComboBox, "StatsChoice")
        StatsChoice.currentIndexChanged.connect(self.stats_chosen)
        StatsChoice.setSeparator("  ")
        StatsChoice.setVisible(False)
        self.WholeOrStat.currentIndexChanged.connect(self.stat_choice)
        self.MultipleRasterCheck.stateChanged.connect(self.stat_choice)
        self.VectorType.currentIndexChanged.connect(self.stat_choice)
        self.Extract.clicked.connect(self.stats_chosen)
        # ADD TO LAYER
        AddToLayer = self.findChild(QtWidgets.QPushButton, 'AddToLayer')
        AddToLayer.clicked.connect(self.add_to_layer)
        ## PROGRESS BAR
        self.extraction_progress.setVisible(False)
        self.extraction_progress.setStyleSheet("QProgressBar::chunk { background-color: #007805; color: #ffe423;}"
                                               "QProgressBar { background-color: #7b7b7b; color: white; text-align: 'center'; \
                                               font-weight: bold;}; border-radius:50px 50px;")

    
    ## ADD RASTER AND OR VECTOR TO LAYER
    def add_to_layer(self):
        file_path =  self.path_to_file(self.PathToFile)
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        folder_path = self.path_to_file(self.PathToFile)

        vector_layer = QgsVectorLayer(vector_path, os.path.splitext(os.path.basename(vector_path))[0], 'ogr')
        raster_layer = QgsRasterLayer(file_path, os.path.splitext(os.path.basename(file_path))[0])
            
        if file_path == "" and vector_path != "" and not self.MultipleRasterCheck.isChecked(): #ok
            QgsProject().instance().addMapLayer(vector_layer)

        elif file_path != "" and vector_path == "" and not self.MultipleRasterCheck.isChecked(): #ok
            QgsProject().instance().addMapLayer(raster_layer)

        elif file_path != "" and vector_path != "" and not self.MultipleRasterCheck.isChecked(): #ok
            QgsProject().instance().addMapLayer(raster_layer)
            QgsProject().instance().addMapLayer(vector_layer)

        elif self.MultipleRasterCheck.isChecked() and os.path.isdir(folder_path) and vector_path == "":
            tiff_files = [os.path.join(folder_path, f) for f in os.listdir(folder_path) if f.endswith(rast_ext)]
            if len(tiff_files) == 0:
                empty_folder()
                return
            raster_layer_folder = QgsRasterLayer(tiff_files[0], os.path.splitext(os.path.basename(tiff_files[0]))[0])
            QgsProject().instance().addMapLayer(raster_layer_folder)

        elif self.MultipleRasterCheck.isChecked() and not os.path.isdir(folder_path) and vector_path != "":
            QgsProject().instance().addMapLayer(vector_layer)

        elif self.MultipleRasterCheck.isChecked() and os.path.isdir(folder_path) and vector_path != "":
            tiff_files = [os.path.join(folder_path, f) for f in os.listdir(folder_path) if f.endswith(rast_ext)]
            if len(tiff_files) == 0:
                empty_folder()
                return
            raster_layer_folder = QgsRasterLayer(tiff_files[0], os.path.splitext(os.path.basename(tiff_files[0]))[0])
            QgsProject().instance().addMapLayer(raster_layer_folder)
            QgsProject().instance().addMapLayer(vector_layer)
        else:
            return
    
    
    def stats_chosen(self):
        choice = self.StatsChoice.currentText().split("  ")
        return choice
    def stat_choice(self):
        if self.MultipleRasterCheck.isChecked() and self.VectorType.currentText() != "Point":
            self.StatsChoice.setVisible(True)
        elif self.SingleRasterCheck.isChecked() and self.WholeOrStat.currentText() == "Zonal statistics" and self.VectorType.currentText() != "Point":
            self.StatsChoice.setVisible(True)
        else:
            self.StatsChoice.setVisible(False)

    # Model adapted to Qt to display data frame from pandas object
    def display_raster_data(self):
        if self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Whole data" and self.SingleRasterCheck.isChecked():
            raster_line = self.raster_line()
            model = PandasModel(raster_line)
            self.raster_table_view.setModel(model)
        elif self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Zonal statistics" and self.SingleRasterCheck.isChecked():
            raster_line_stats = self.raster_line_stats()
            model = PandasModel(raster_line_stats)
            self.raster_table_view.setModel(model)
        elif self.VectorType.currentText() == "Line" and self.MultipleRasterCheck.isChecked():
            multi_raster_line_stats = self.multi_raster_line_stats()
            model = PandasModel(multi_raster_line_stats)
            self.raster_table_view.setModel(model)
        elif self.VectorType.currentText() == "Point" and self.SingleRasterCheck.isChecked(): #self.SingleRasterCheck.isChecked():
            with_point = self.raster_point_df()
            model = PandasModel(with_point)
            self.raster_table_view.setModel(model)
            ### MULTI RASTER
        elif self.VectorType.currentText() == "Point" and self.MultipleRasterCheck.isChecked(): #self.SingleRasterCheck.isChecked():
            multi_raster_point = self.multi_raster_point()
            model = PandasModel(multi_raster_point)
            self.raster_table_view.setModel(model)
        elif self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Whole data" and not self.MultipleRasterCheck.isChecked():
            whole_zonal_value = self.whole_polygon()
            model = PandasModel(whole_zonal_value)
            self.raster_table_view.setModel(model)
        elif self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Zonal statistics" and self.SingleRasterCheck.isChecked():
            single_zonal_stat = self.single_zonal_stat()
            model = PandasModel(single_zonal_stat)
            self.raster_table_view.setModel(model)
            ### MULTIPLE RASTER -- FOLDER
        elif self.VectorType.currentText() == "Polygon" and self.MultipleRasterCheck.isChecked():
            multi_zonal_stat = self.mult_zonal_stat()
            model = PandasModel(multi_zonal_stat)
            self.raster_table_view.setModel(model)
        elif self.SingleRasterCheck.isChecked() and self.VectorType.currentIndex() not in [1,2,3]:
            sing_raster = self.single_raster_df()
            model = PandasModel(sing_raster)
            self.raster_table_view.setModel(model)
        self.file_path = None

    # SINGLE RASTER
    def single_raster_df(self):
        file_path = self.path_to_file(self.PathToFile)
        if not file_path and self.SingleRasterCheck.isChecked() and self.VectorType.currentIndex() == 0:
            missed_raster_path()
            return
        elif os.path.isdir(file_path):
            is_directory()
            return
        else:
            ds = gdal.Open(file_path)
            # Get the transform from pixel coordinates to georeferenced coordinates
            transform = ds.GetGeoTransform()
            # Get the number of rows and columns in the raster
            rows = ds.RasterYSize
            cols = ds.RasterXSize
            # Get the first band of the raster
            band = ds.GetRasterBand(1)
            # Create an empty list to store the longitudes, latitudes, and values
            lons_lats_vals = []
            # Loop over all the pixels in the raster
            for row in range(rows):
                for col in range(cols):
                    # Get the longitude and latitude of the current pixel
                    lon = transform[0] + col * transform[1] + row * transform[2]
                    lat = transform[3] + col * transform[4] + row * transform[5]      
                    # Read the value of the current pixel
                    value = np.round(band.ReadAsArray(col, row, 1, 1)[0][0], 3)
                    # Add the longitude, latitude, and value to the list
                    lons_lats_vals.append((lon, lat, value))
            # Convert the list to a pandas data frame
            RasterDf = pd.DataFrame(lons_lats_vals, columns=["longitude", "latitude", "value"])
            RasterDf = RasterDf[~np.isinf(RasterDf['value'])]
            RasterDf['value'] = RasterDf['value']
            df = RasterDf.reset_index()
            del df["index"]
            df.index += 1
            return df
    
    # RASTER AND POINT
    def raster_point_df(self):
        file_path = self.path_to_file(self.PathToFile)
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        if not vector_path and not file_path and self.VectorType.currentText() == "Point":
            missed_raster_vector()
            return
        elif not vector_path and self.VectorType.currentText() == "Point":
            missed_vector_path()
            return
        elif not file_path and self.VectorType.currentText() == "Point":
            missed_raster_path()
            return
        elif os.path.isdir(file_path):
            is_directory()
            return
        elif self.handle_shp_feature_type_error()  != "POINT" and self.handle_shp_feature_type_error()  != "MULTIPOINT":
            not_point()
            return
        else:
            # Open the vector layer
            vector_ds = ogr.Open(vector_path)
            layer = vector_ds.GetLayer()
            vector_sr = layer.GetSpatialRef()
            # Load the raster layer
            raster_ds = gdal.Open(file_path)
            raster_band = raster_ds.GetRasterBand(1)
            # Get the CRS of the raster layer
            raster_sr = osr.SpatialReference()
            raster_sr.ImportFromWkt(raster_ds.GetProjection())
            values = []
            ## create transformer
            transformer = osr.CoordinateTransformation(vector_sr, raster_sr)
            if "None" in str(transformer):
                Unoverlable()
                return
            else:
                try:
                    for feature in layer:
                        geometry = feature.GetGeometryRef()
                        geometry.Transform(transformer)

                        if geometry.GetGeometryName() == 'MULTIPOINT':
                            for i in range(geometry.GetGeometryCount()):
                                point = geometry.GetGeometryRef(i)
                                x, y = point.GetX(), point.GetY()
                                # Get the raster value for the point
                                transform = raster_ds.GetGeoTransform()
                                px = int((x - transform[0]) / transform[1])
                                py = int((y - transform[3]) / transform[5])
                                raster_band = raster_ds.GetRasterBand(1)
                                raster_array = raster_band.ReadAsArray(px, py, 1, 1)[0,0]
                                values.append((feature.GetField(self.field_chosen()), x, y, 
                                            round(float(raster_array), 7)))
                        else:
                            x, y = geometry.GetX(), geometry.GetY()
                            # Get the raster value for the point
                            transform = raster_ds.GetGeoTransform()
                            px = int((x - transform[0]) / transform[1])
                            py = int((y - transform[3]) / transform[5])
                            raster_band = raster_ds.GetRasterBand(1)
                            raster_array = raster_band.ReadAsArray(px, py, 1, 1)[0,0]
                            values.append((feature.GetField(self.field_chosen()), x, y, 
                                        round(float(raster_array), 3)))
                except:
                    pass
                df = pd.DataFrame(values, columns = [self.field_chosen(), "longitude", "latitude", "value"])
                df.index += 1
                return df
    
    # MULTI RASTER AND POINT
    def multi_raster_point(self):
        folder_path = self.path_to_file(self.PathToFile)
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        try:
            tiff_files = [os.path.join(folder_path, f) for f in os.listdir(folder_path) if f.endswith(rast_ext)]
        except:
            pass
        values = []
        if not vector_path and self.VectorType.currentText() == "Point" and self.MultipleRasterCheck.isChecked():
            missed_vector_path()
            return
        elif self.MultipleRasterCheck.isChecked() and not os.path.isdir(folder_path):
            no_folder()
            return
        elif len(tiff_files) == 0:
            empty_folder()
            return
        elif self.handle_shp_feature_type_error()  != "POINT" and self.handle_shp_feature_type_error()  != "MULTIPOINT":
            not_point()
            return
        else:
            # Open the vector layer
            vector_ds = ogr.Open(vector_path)
            layer = vector_ds.GetLayer()
            vector_srs = layer.GetSpatialRef()
            # set progress bar
            self.extraction_progress.setVisible(True)
            self.extraction_progress.setMinimum(0) 
            self.extraction_progress.setMaximum(len(tiff_files))  # set maximum value to number of raster files
            for i, file in enumerate(tiff_files):
                self.extraction_progress.setValue(i+1)  # increment progress bar value
                # get name of raster file without extension
                name = os.path.splitext(os.path.basename(file))[0]
                # Open the raster file
                raster = gdal.Open(file)
                geo_transform = raster.GetGeoTransform()
                proj = raster.GetProjection()
                raster_src = osr.SpatialReference()
                raster_src.ImportFromWkt(proj)
                transform = osr.CoordinateTransformation(vector_srs, raster_src)
                if "None" in str(transform):
                    Unoverlable()
                    return
                # loop through tif files in folder
                for feature in layer:
                    ID = feature.GetField(self.field_chosen())
                    geom = feature.GetGeometryRef()
                    geom.Transform(transform)

                    try:
                        if geom.GetGeometryName() == 'MULTIPOINT':
                            for i in range(geom.GetGeometryCount()):
                                point = geom.GetGeometryRef(i)
                                x, y = point.GetX(), point.GetY()
                                px = int((x - geo_transform[0]) / geo_transform[1])
                                py = int((y - geo_transform[3]) / geo_transform[5])
                                value = raster.GetRasterBand(1).ReadAsArray(px, py, 1, 1)[0][0]
                                values.append({"source":name, self.field_chosen():ID, "values": value})
                        else:
                            point = np.array(geom.GetPoints())
                            for i in range(len(point)):
                                x, y = point[i]
                                px = int((x - geo_transform[0]) / geo_transform[1])
                                py = int((y - geo_transform[3]) / geo_transform[5])
                                value = raster.GetRasterBand(1).ReadAsArray(px, py, 1, 1)[0][0]
                                values.append({"source":name, self.field_chosen():ID, "values": value})
                    except:
                        pass

            self.extraction_progress.setValue(len(tiff_files))  # set progress bar to full
            # Create a Pandas DataFrame to store the extracted values
            df = pd.DataFrame(values)
            df.index += 1
            global multi_raster_point
            multi_raster_point = df 
            return multi_raster_point

    # SINGLE ZONAL COMPUTATION
    def single_zonal_stat(self):
        file_path = self.path_to_file(self.PathToFile)
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        stat_name = self.stats_chosen()
        if not vector_path and not file_path and self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Zonal statistics":
            missed_raster_vector()
            return
        elif not vector_path and self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Zonal statistics":
            missed_vector_path()
            return
        elif not file_path and self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Zonal statistics":
            missed_raster_path()
            return
        elif os.path.isdir(file_path):
            is_directory()
            return
        elif self.handle_shp_feature_type_error() != "POLYGON" and self.handle_shp_feature_type_error() != "MULTIPOLYGON":
            not_polygon()
            return
        elif stat_name == [""]:
            no_stat()
            return

        else:
            # Open the vector layer
            vector_ds = ogr.Open(vector_path)
            layer = vector_ds.GetLayer()
            vector_sr = layer.GetSpatialRef()
            # Load the raster layer
            raster_ds = gdal.Open(file_path)
            raster_sr = osr.SpatialReference()
            raster_sr.ImportFromWkt(raster_ds.GetProjection())
            ## create transformer
            transformer = osr.CoordinateTransformation(vector_sr, raster_sr)
            print(transformer)
            if "None" in str(transformer):
                Unoverlable()
                return
            out_df = zonal_data(file_path, vector_path, self.field_chosen())
            agg_df = groupby_agg(df=out_df, group_cols = self.field_chosen(), value_col ="value", funcs = stat_name)
            agg_df.rename(columns={'{}'.format(self.field_chosen()): self.field_chosen()}, inplace=True)
            return agg_df
    
    # RASTER AND LINE
    def raster_line(self):
        file_path = self.path_to_file(self.PathToFile)
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        if not vector_path and not file_path and self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Whole data":
            missed_raster_vector()
            return
        elif not vector_path and self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Whole data":
            missed_vector_path()
            return
        elif not file_path and self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Whole data":
            missed_raster_path()
            return
        elif os.path.isdir(file_path):
            is_directory()
            return
        elif self.handle_shp_feature_type_error() != "LINESTRING" and self.WholeOrStat.currentText() == "Whole data":
            not_line()
            return
        else:
            # Open the raster file
            raster = gdal.Open(file_path)
            geo_transform = raster.GetGeoTransform()
            proj = raster.GetProjection()
            raster_src = osr.SpatialReference()
            raster_src.ImportFromWkt(proj)

            # Open the vector file
            vector = ogr.Open(vector_path)
            layer = vector.GetLayer()
            vector_srs = layer.GetSpatialRef()

            transform = osr.CoordinateTransformation(vector_srs, raster_src)
            if "None" in str(transform):
                Unoverlable()
                return
            # Initialize an empty list to store the extracted values
            values = []
            # Loop through the features in the input vector layer and extract the values along each feature
            try:
                for feature in layer:
                    geom = feature.GetGeometryRef()
                    geom.Transform(transform)
                    line = np.array(geom.GetPoints())
                    #print(line)
                    for i in range(len(line)):
                        x, y = line[i]
                        px = int((x - geo_transform[0]) / geo_transform[1])
                        py = int((y - geo_transform[3]) / geo_transform[5])
                        value = raster.GetRasterBand(1).ReadAsArray(px, py, 1, 1)[0][0]
                        values.append({self.field_chosen():feature.GetField(self.field_chosen()),
                                    "longitude":x, "latitude":y, "values": round(float(value), 3)})
            except:
                pass
            df = pd.DataFrame(values)
            df.index += 1
            return df

    # RASTER AND LINE STATS
    def raster_line_stats(self):
        file_path = self.path_to_file(self.PathToFile)
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        stat_name = self.stats_chosen()
        if not vector_path and not file_path and self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Zonal statistics":
            missed_raster_vector()
            return
        elif not vector_path and self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Zonal statistics":
            missed_vector_path()
            return
        elif not file_path and self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Zonal statistics":
            missed_raster_path()
            return
        elif os.path.isdir(file_path):
            is_directory()
            return
        elif self.handle_shp_feature_type_error() != "LINESTRING" and self.WholeOrStat.currentText() == "Zonal statistics":
            not_line()
            return
        elif stat_name == [""]:
            no_stat()
            return
        else:
            # Open the raster file
            raster = gdal.Open(file_path)
            geo_transform = raster.GetGeoTransform()
            proj = raster.GetProjection()
            raster_src = osr.SpatialReference()
            raster_src.ImportFromWkt(proj)

            # Open the vector file
            vector = ogr.Open(vector_path)
            layer = vector.GetLayer()
            vector_srs = layer.GetSpatialRef()

            transform = osr.CoordinateTransformation(vector_srs, raster_src)
            if "None" in str(transform):
                Unoverlable()
                return
            # Initialize an empty list to store the extracted values
            values = []
            stat_name = self.stats_chosen()
            # Loop through the features in the input vector layer and extract the values along each feature
            try:
                for feature in layer:
                    geom = feature.GetGeometryRef()
                    geom.Transform(transform)
                    line = np.array(geom.GetPoints())
                    #print(line)
                    for i in range(len(line)):
                        x, y = line[i]
                        px = int((x - geo_transform[0]) / geo_transform[1])
                        py = int((y - geo_transform[3]) / geo_transform[5])
                        value = raster.GetRasterBand(1).ReadAsArray(px, py, 1, 1)[0][0]
                        values.append({self.field_chosen():feature.GetField(self.field_chosen()),
                                    "val": float(value)})
            except:
                pass
            result_df = pd.DataFrame(values)
            df = groupby_agg(df = result_df, group_cols = self.field_chosen(), 
                             value_col = "val", funcs = stat_name)
            return df

   # MULTI RASTER AND LINE STATS
    def multi_raster_line_stats(self):
        folder_path = self.path_to_file(self.PathToFile)
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        stat_name = self.stats_chosen()
        try:
            tiff_files = [os.path.join(folder_path, f) for f in os.listdir(folder_path) if f.endswith(rast_ext)]
        except:
            pass
        values = []
        if not vector_path and self.VectorType.currentText() == "Line" and self.MultipleRasterCheck.isChecked():
            missed_vector_path()
            return
        elif self.MultipleRasterCheck.isChecked() and not os.path.isdir(folder_path):
            no_folder()
            return
        elif len(tiff_files) == 0:
            empty_folder()
            return
        elif self.handle_shp_feature_type_error()  != "LINESTRING":
            not_line()
            return
        elif stat_name == [""]:
            no_stat()
            return
        else:
            # Open the raster file
            shp = ogr.Open(vector_path)
            layer = shp.GetLayer()
            vector_srs = layer.GetSpatialRef()
            values = []
            stat_name = self.stats_chosen()
            ## set progress bar
            self.extraction_progress.setVisible(True)
            self.extraction_progress.setMinimum(0) 
            self.extraction_progress.setMaximum(len(tiff_files))  # set maximum value to number of raster files
            for i, file in enumerate(tiff_files):
                self.extraction_progress.setValue(i+1)  # increment progress bar value
                # get name of raster file without extension
                name = os.path.splitext(os.path.basename(file))[0]
                # Open the raster file
                raster = gdal.Open(file)
                geo_transform = raster.GetGeoTransform()
                proj = raster.GetProjection()
                raster_src = osr.SpatialReference()
                raster_src.ImportFromWkt(proj)
                transform = osr.CoordinateTransformation(vector_srs, raster_src)
                if "None" in str(transform):
                    Unoverlable()
                    return
            # loop through tif files in folder
                try:
                    for feature in layer:
                        ID = feature.GetField(self.field_chosen())
                        geom = feature.GetGeometryRef()
                        geom.Transform(transform)
                        line = np.array(geom.GetPoints())
                        for i in range(len(line)):
                            x, y = line[i]
                            px = int((x - geo_transform[0]) / geo_transform[1])
                            py = int((y - geo_transform[3]) / geo_transform[5])
                            value = raster.GetRasterBand(1).ReadAsArray(px, py, 1, 1)[0][0]
                            values.append({"source":name, self.field_chosen():ID, "values": value})
                except:
                    pass
            self.extraction_progress.setValue(len(tiff_files))  # set progress bar to full
            # Create a Pandas DataFrame to store the extracted values
            result_df = pd.DataFrame(values)
            df = groupby_agg(df = result_df, group_cols = ["source", self.field_chosen()], 
                             value_col = "values", funcs = stat_name)
            global multi_raster_line_stats
            multi_raster_line_stats = df
            return multi_raster_line_stats
 
    # WHOLE DATA
    def whole_polygon(self):
        file_path = self.path_to_file(self.PathToFile)
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        if not vector_path and not file_path and self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Whole data":
            missed_raster_vector()
            return
        elif not vector_path and self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Whole data":
            missed_vector_path()
            return
        elif not file_path and self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Whole data":
            missed_raster_path()
            return
        elif os.path.isdir(file_path):
            is_directory()
            return
        elif self.handle_shp_feature_type_error() != "POLYGON" and self.handle_shp_feature_type_error() != "MULTIPOLYGON" and self.WholeOrStat.currentText() == "Whole data":
            not_polygon()
            return
        else:
            out_df = zonal_data(file_path, vector_path, self.field_chosen())
            return out_df

    # MULTIPLE ZONAL STATS
    def mult_zonal_stat(self):
        self.extraction_progress.setVisible(True)
        folder_path = self.path_to_file(self.PathToFile)
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        stat_name = self.stats_chosen()
        # create empty dataframe for output
        output_list = []
        try:
            tiff_files = [os.path.join(folder_path, f) for f in os.listdir(folder_path) if f.endswith(rast_ext)]
        except:
            pass
        if not vector_path and self.VectorType.currentText() == "Polygon" and self.MultipleRasterCheck.isChecked():
            missed_vector_path()
            return

        elif not os.path.isdir(folder_path) and self.VectorType.currentText() == "Polygon" and self.MultipleRasterCheck.isChecked():
            no_folder()
            return
        elif len(tiff_files) == 0:
            empty_folder()
            return
        elif self.handle_shp_feature_type_error() != "POLYGON" and self.handle_shp_feature_type_error() != "MULTIPOLYGON":
            not_polygon()
            return
        elif len(tiff_files) == 0:
            empty_folder()
            return
        elif stat_name == [""]:
            no_stat()
            return
        else:
            self.extraction_progress.setVisible(True)
            self.extraction_progress.setMinimum(0) 
            self.extraction_progress.setMaximum(len(tiff_files)) 
            for  i, file in enumerate(tiff_files):
                self.extraction_progress.setValue(i+1)  # increment progress bar value
                # get name of raster file without extension
                source = os.path.splitext(os.path.basename(file))[0]
                out_df = zonal_data(file, vector_path, self.field_chosen())
                agg_df = groupby_agg(df=out_df, group_cols=self.field_chosen(), value_col="value", funcs = stat_name)
                agg_df["source"] = source
                output_list.append(agg_df)
            self.extraction_progress.setValue(len(tiff_files))  # set progress bar to full
            output = pd.concat(output_list, ignore_index=True)
            output.rename(columns={'{}'.format(self.field_chosen()): self.field_chosen()}, inplace=True)
            output.insert(0, "source", output.pop("source"))
            global mult_zonal_stat
            mult_zonal_stat = output
            return mult_zonal_stat
            
    
    # CLEAR TABLE FUNCTION 
    def clear_raster_table(self):
        self.extraction_progress.setVisible(False)
        model = self.raster_table_view.model()
        model.removeRows(0, model.rowCount())
        self.raster_table_view.setModel(PandasModel(self.originalDataFrame))

    # SAVE TABLE
    def save_raster_table(self):
        # Show a file dialog to choose the file to save the data to
        fileDialog = QtWidgets.QFileDialog()
        fileDialog.setAcceptMode(QtWidgets.QFileDialog.AcceptSave)
        fileDialog.setDefaultSuffix('csv')
        fileDialog.setNameFilter('CSV Files (*.csv);; Text (*.txt)')

        # save with point
        if self.VectorType.currentText() == "Point" and self.SingleRasterCheck.isChecked():
            if fileDialog.exec_():
                fileName = fileDialog.selectedFiles()[0]
                self.raster_point_df().to_csv(fileName, index=False)
        # MULTIPLE RASTER -- POINT
        if self.VectorType.currentText() == "Point" and self.MultipleRasterCheck.isChecked():
            if fileDialog.exec_():
                fileName = fileDialog.selectedFiles()[0]
                multi_raster_point.to_csv(fileName, index=False)
        # save zingle zonal
        elif self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Zonal statistics" and self.SingleRasterCheck.isChecked():
            if fileDialog.exec_():
                fileName = fileDialog.selectedFiles()[0]
                self.single_zonal_stat().to_csv(fileName, index=False)
        # save whole zonal value
        elif self.VectorType.currentText() == "Polygon" and self.WholeOrStat.currentText() == "Whole data" and self.SingleRasterCheck.isChecked():
            if fileDialog.exec_():
                fileName = fileDialog.selectedFiles()[0]
                self.whole_polygon().to_csv(fileName, index=False)
        ### MULTIPLE RASTER -- FOLDER
        elif self.VectorType.currentText() == "Polygon" and self.MultipleRasterCheck.isChecked():
            if fileDialog.exec_():
                fileName = fileDialog.selectedFiles()[0]
                mult_zonal_stat.to_csv(fileName, index=False)
        # save line
        elif self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Whole data" and self.SingleRasterCheck.isChecked():
            if fileDialog.exec_():
                fileName = fileDialog.selectedFiles()[0]
                self.raster_line().to_csv(fileName, index=False)
        elif self.VectorType.currentText() == "Line" and self.WholeOrStat.currentText() == "Zonal statistics" and self.SingleRasterCheck.isChecked():
            if fileDialog.exec_():
                fileName = fileDialog.selectedFiles()[0]
                self.raster_line_stats().to_csv(fileName, index=False)
        elif self.VectorType.currentText() == "Line" and self.MultipleRasterCheck.isChecked():
            if fileDialog.exec_():
                fileName = fileDialog.selectedFiles()[0]
                multi_raster_line_stats.to_csv(fileName, index=False)
        elif self.SingleRasterCheck.isChecked() and self.VectorType.currentIndex() not in [1,2,3]:
            if fileDialog.exec_():
                fileName = fileDialog.selectedFiles()[0]
                self.single_raster_df().to_csv(fileName, index=False)

    # ENABLE SAVE BUTTON
    def enableSave(self):
        mod = self.raster_table_view.model()
        if mod.rowCount() != 0:
            self.Save.setEnabled(True)
            self.Clear.setEnabled(True)

    def unableSave(self):
        mod = self.raster_table_view.model()
        self.Save.setEnabled(False)
        if mod.rowCount() != 0:
            self.Save.setEnabled(False)


 
    def path_to_file(self, PathToFile):
        if self.SingleRasterCheck.isChecked():
            self.PathToFile.setStorageMode(QgsFileWidget.GetFile)
            path = self.PathToFile.filePath()
            return path
        elif self.MultipleRasterCheck.isChecked():
            self.PathToFile.setStorageMode(QgsFileWidget.GetDirectory)
            path = self.PathToFile.filePath()
            return path

    ### VECTORE FILE PATH
    def vector_file_path(self, PathToVectorFile):
            #global vector_path
            vector_path = PathToVectorFile.filePath()
            
            if vector_path.endswith("csv") and self.VectorType.currentText() == "Point":
               # To use the dialog, instantiate and show it
               try:
                   dialog = Delim(file_path=vector_path)
                   dialog.exec_()
                   a, y = os.path.splitext(os.path.basename(vector_path))
                   vector_path = vector_path.replace(y, ".shp")
                   PathToVectorFile.lineEdit().setText(vector_path)
               except:
                   pass

            return vector_path
    
    # DISCONNECT ALL FILE/FOLDER CHOSED
    def remove_whole_or_stat(self):
        self.WholeOrStat.setVisible(False)

    ##  VECTOR FILE CHOICE APPEARING 
    def selectVectorFile(self, index):
        if not index == 0:
            self.PathToVectorFile.setVisible(True)
            self.FieldName.setVisible(True)
        #
        else:
            self.PathToVectorFile.setVisible(False)
            self.FieldName.setVisible(False)
    
    # WHOLE OR STAT APPEARING
    def whole_or_stat(self, index):
        if index in [1,3] and self.SingleRasterCheck.isChecked():
            self.WholeOrStat.setVisible(True)
        else:
            self.WholeOrStat.setVisible(False)

    ## VECTOR FIELD NAME
    def field_name(self):
        self.FieldName.clear()
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        try:
            vector_ds = ogr.Open(vector_path)
            layer = vector_ds.GetLayer()
            layer_defn = layer.GetLayerDefn()
            shp_colname = []
            # Loop through the field definitions and get each field
            for i in range(layer_defn.GetFieldCount()):
                layer_defn = layer.GetLayerDefn()
                names = layer_defn.GetFieldDefn(i).GetName()
                shp_colname.append(names)
            self.FieldName.addItems(shp_colname)
        except:
            pass
    
    
    ## VECTOR FIELD CHOSEN
    def field_chosen(self):
        field_chosen = self.FieldName.currentText()
        return field_chosen
    
    ## HANDLE SHP FEATURE TYPE ERROR
    def handle_shp_feature_type_error(self):
        vector_path = r"{}".format(self.vector_file_path(self.PathToVectorFile))
        try:
            dataset = ogr.Open(vector_path)
            layer = dataset.GetLayer()
            for feature in layer:
                geometry = feature.GetGeometryRef()
                feature_name = geometry.GetGeometryName()
            return str(feature_name)
        except:
            pass
    
            
#####################