Fix bug/194 - Raster identification tool

ee_plugin/provider.py

@@ -5,6 +5,7 @@
 
 import json
 
+import ee
 from qgis.core import (
     Qgis,
     QgsCoordinateReferenceSystem,
@@ -20,6 +21,8 @@
 )
 from qgis.PyQt.QtCore import QObject
 
+from ee_plugin import Map
+
 BAND_TYPES = {
     "int8": Qgis.Int16,
     "int16": Qgis.Int16,
@@ -239,18 +242,14 @@ def htmlMetadata(self):
     def identify(
         self, point, format, boundingBox=None, width=None, height=None, dpi=None
     ):
-        # TODO: speed-up, extend this to maintain cache of visible image, update cache on-the-fly when needed
-        import ee
-
-        from ee_plugin import Map, utils
+        dataset_projection = self.ee_info["bands"][0]["crs"]
+        point_ee = ee.Geometry.Point(
+            [point.x(), point.y()], self.crs().authid()
+        ).transform(dataset_projection)
 
-        point = utils.geom_to_geo(point)
-        point_ee = ee.Geometry.Point([point.x(), point.y()])
+        reducer = ee.Reducer.first()
         scale = Map.getScale()
-        value = self.ee_object.reduceRegion(
-            ee.Reducer.first(), point_ee, scale
-        ).getInfo()
-
+        value = self.ee_object.reduceRegion(reducer, point_ee, scale).getInfo()
         band_indices = range(1, self.bandCount() + 1)
         band_names = [self.generateBandName(band_no) for band_no in band_indices]
         band_values = [value[band_name] for band_name in band_names]
@@ -260,8 +259,6 @@ def identify(
 
         return result
 
-    # sample()
-
     def lastErrorTitle(self):
         return self.wms.lastErrorTitle()
 
@@ -427,17 +424,14 @@ def set_ee_object(self, ee_object):
 
 
 class EarthEngineVectorDataProvider(QgsVectorDataProvider):
-    # TODO
     pass
 
 
 class EarthEngineRasterCollectionDataProvider(QgsRasterDataProvider):
-    # TODO
     pass
 
 
 class EarthEngineVectorCollectionDataProvider(QgsVectorDataProvider):
-    # TODO
     pass
 
 

ee_plugin/utils.py

@@ -8,12 +8,8 @@
 import ee
 import qgis
 from qgis.core import (
-    QgsCoordinateReferenceSystem,
-    QgsCoordinateTransform,
-    QgsPointXY,
     QgsProject,
     QgsRasterLayer,
-    QgsRectangle,
 )
 
 from .ee_plugin import VERSION as ee_plugin_version
@@ -183,16 +179,3 @@ def add_ee_catalog_image(name, asset_name, visParams, collection_props):
 def check_version():
     # check if we have the latest version only once plugin is used, not once it is loaded
     qgis.utils.plugins["ee_plugin"].check_version()
-
-
-def geom_to_geo(geom):
-    crs_src = QgsCoordinateReferenceSystem(QgsProject.instance().crs())
-    crs_dst = QgsCoordinateReferenceSystem("EPSG:4326")
-    proj2geo = QgsCoordinateTransform(crs_src, crs_dst, QgsProject.instance())
-
-    if isinstance(geom, QgsPointXY):
-        return proj2geo.transform(geom)
-    elif isinstance(geom, QgsRectangle):
-        return proj2geo.transformBoundingBox(geom)
-    else:
-        return geom.transform(proj2geo)

test/test_provider.py

@@ -0,0 +1,65 @@
+from unittest.mock import patch
+
+import ee
+from qgis.core import QgsProject, QgsPointXY
+
+from ee_plugin import Map
+
+
+def test_raster_identify():
+    # Add an image layer to QGIS
+    image = ee.Image("USGS/SRTMGL1_003")
+    vis_params = {
+        "min": 0,
+        "max": 4000,
+        "palette": ["006633", "E5FFCC", "662A00", "D8D8D8", "F5F5F5"],
+    }
+
+    # mocking Map.getScale() was necessary for ee to reproject point correctly
+    # simpler than configuring iface and other QGIS objects
+    mock_get_scale = patch("ee_plugin.Map.getScale", return_value=625.7583)
+    mock_get_scale.start()
+
+    Map.addLayer(image, vis_params, "DEM")
+
+    qgis_instance = QgsProject.instance()
+    layer = qgis_instance.mapLayersByName("DEM")[0]
+    assert layer is not None, "Layer 'DEM' was not added to the project"
+    assert layer.crs().authid() == "EPSG:3857", "Layer CRS does not match project CRS"
+
+    provider = layer.dataProvider()
+    qgis_point = QgsPointXY(-13551778.88787266425788403, 5917193.28679858986288309)
+    raster_identify_result = provider.identify(
+        qgis_point, format=1, height=1, width=1, dpi=96
+    )
+
+    assert raster_identify_result, "Identify operation returned no results"
+    assert (
+        raster_identify_result.results()
+    ), "Identify operation returned an empty result set"
+    assert (
+        raster_identify_result.results()[1] > 0
+    ), "Identified elevation is not positive"
+
+
+def test_reduce_region():
+    ee.Initialize()
+
+    image = ee.Image("USGS/SRTMGL1_003")
+    vis_params = {
+        "min": 0,
+        "max": 4000,
+        "palette": ["006633", "E5FFCC", "662A00", "D8D8D8", "F5F5F5"],
+    }
+    Map.addLayer(image, vis_params, "DEM")
+    Map.setCenter(-121.753, 46.855, 9)
+
+    test_point = ee.Geometry.Point([-121.753, 46.855], "EPSG:4326")
+    scale = 30  # Scale in meters
+    reducer = ee.Reducer.first()
+    result = image.reduceRegion(
+        reducer=reducer, geometry=test_point, scale=scale
+    ).getInfo()
+
+    assert "elevation" in result, {"message": "Elevation not found in result"}
+    assert result["elevation"] > 0, {"message": "Elevation is not positive"}