EmileSonneveld
commented
7 months ago Graph to test the issue:
{
"process_graph": {
"loadcollection1": {
"process_id": "load_collection",
"arguments": {
"bands": [
"B02",
"B03",
"B04",
"B05",
"B06",
"B07",
"B08",
"B11",
"B12",
"SCL"
],
"id": "SENTINEL2_L2A",
"properties": {
"tileId": {
"process_graph": {
"eq1": {
"process_id": "eq",
"arguments": {
"x": {
"from_parameter": "value"
},
"y": "33UV*"
},
"result": true
}
}
},
"eo:cloud_cover": {
"process_graph": {
"lte1": {
"process_id": "lte",
"arguments": {
"x": {
"from_parameter": "value"
},
"y": 95
},
"result": true
}
}
}
},
"spatial_extent": {
"east": 4660000,
"south": 3240000,
"west": 4640000,
"north": 3260000,
"crs": "EPSG:3035"
},
"temporal_extent": [
"2018-03-01",
"2018-10-31"
]
}
},
"resamplespatial1": {
"process_id": "resample_spatial",
"arguments": {
"align": "upper-left",
"data": {
"from_node": "loadcollection1"
},
"method": "near",
"projection": 3035,
"resolution": 10
}
},
"maskscldilation1": {
"process_id": "mask_scl_dilation",
"arguments": {
"data": {
"from_node": "resamplespatial1"
},
"erosion_kernel_size": 3,
"kernel1_size": 17,
"kernel2_size": 77,
"mask1_values": [
2,
4,
5,
6,
7
],
"mask2_values": [
3,
8,
9,
10,
11
],
"scl_band_name": "SCL"
}
},
"filterbands1": {
"process_id": "filter_bands",
"arguments": {
"bands": [
"B02",
"B03",
"B04",
"B05",
"B06",
"B07",
"B08",
"B11",
"B12"
],
"data": {
"from_node": "maskscldilation1"
}
}
},
"aggregatetemporalperiod1": {
"process_id": "aggregate_temporal_period",
"arguments": {
"data": {
"from_node": "filterbands1"
},
"period": "dekad",
"reducer": {
"process_graph": {
"median1": {
"process_id": "median",
"arguments": {
"data": {
"from_parameter": "data"
}
},
"result": true
}
}
}
}
},
"applydimension1": {
"process_id": "apply_dimension",
"arguments": {
"data": {
"from_node": "aggregatetemporalperiod1"
},
"dimension": "t",
"process": {
"process_graph": {
"arrayinterpolatelinear1": {
"process_id": "array_interpolate_linear",
"arguments": {
"data": {
"from_parameter": "data"
}
},
"result": true
}
}
}
}
},
"apply1": {
"process_id": "apply",
"arguments": {
"data": {
"from_node": "applydimension1"
},
"process": {
"process_graph": {
"linearscalerange1": {
"process_id": "linear_scale_range",
"arguments": {
"inputMax": 65534,
"inputMin": 0,
"outputMax": 65534,
"outputMin": 0,
"x": {
"from_parameter": "x"
}
},
"result": true
}
}
}
}
},
"loadcollection2": {
"process_id": "load_collection",
"arguments": {
"bands": [
"VH",
"VV"
],
"id": "SENTINEL1_GRD",
"properties": {
"sat:orbit_state": {
"process_graph": {
"eq2": {
"process_id": "eq",
"arguments": {
"x": {
"from_parameter": "value"
},
"y": "DESCENDING"
},
"result": true
}
}
},
"polarisation": {
"process_graph": {
"eq3": {
"process_id": "eq",
"arguments": {
"x": {
"from_parameter": "value"
},
"y": "VV&VH"
},
"result": true
}
}
}
},
"spatial_extent": {
"east": 4660000,
"south": 3240000,
"west": 4640000,
"north": 3260000,
"crs": "EPSG:3035"
},
"temporal_extent": [
"2018-03-01",
"2018-10-31"
]
}
},
"sarbackscatter1": {
"process_id": "sar_backscatter",
"arguments": {
"coefficient": "sigma0-ellipsoid",
"contributing_area": false,
"data": {
"from_node": "loadcollection2"
},
"elevation_model": "COPERNICUS_30",
"ellipsoid_incidence_angle": false,
"local_incidence_angle": false,
"mask": false,
"noise_removal": true,
"options": {
"implementation_version": "2",
"tile_size": 256,
"otb_memory": 1024,
"debug": false,
"elev_geoid": "/opt/openeo-vito-aux-data/egm96.tif"
}
}
},
"resamplespatial2": {
"process_id": "resample_spatial",
"arguments": {
"align": "upper-left",
"data": {
"from_node": "sarbackscatter1"
},
"method": "near",
"projection": 3035,
"resolution": 10
}
},
"aggregatetemporalperiod2": {
"process_id": "aggregate_temporal_period",
"arguments": {
"data": {
"from_node": "resamplespatial2"
},
"period": "dekad",
"reducer": {
"process_graph": {
"mean1": {
"process_id": "mean",
"arguments": {
"data": {
"from_parameter": "data"
}
},
"result": true
}
}
}
}
},
"applydimension2": {
"process_id": "apply_dimension",
"arguments": {
"data": {
"from_node": "aggregatetemporalperiod2"
},
"dimension": "t",
"process": {
"process_graph": {
"arrayinterpolatelinear2": {
"process_id": "array_interpolate_linear",
"arguments": {
"data": {
"from_parameter": "data"
}
},
"result": true
}
}
}
}
},
"applydimension3": {
"process_id": "apply_dimension",
"arguments": {
"data": {
"from_node": "applydimension2"
},
"dimension": "bands",
"process": {
"process_graph": {
"arrayelement1": {
"process_id": "array_element",
"arguments": {
"data": {
"from_parameter": "data"
},
"index": 0
}
},
"log1": {
"process_id": "log",
"arguments": {
"base": 10,
"x": {
"from_node": "arrayelement1"
}
}
},
"multiply1": {
"process_id": "multiply",
"arguments": {
"x": 10,
"y": {
"from_node": "log1"
}
}
},
"add1": {
"process_id": "add",
"arguments": {
"x": {
"from_node": "multiply1"
},
"y": 83
}
},
"divide1": {
"process_id": "divide",
"arguments": {
"x": {
"from_node": "add1"
},
"y": 20
}
},
"power1": {
"process_id": "power",
"arguments": {
"base": 10,
"p": {
"from_node": "divide1"
}
}
},
"arrayelement2": {
"process_id": "array_element",
"arguments": {
"data": {
"from_parameter": "data"
},
"index": 0
}
},
"isnodata1": {
"process_id": "is_nodata",
"arguments": {
"x": {
"from_node": "arrayelement2"
}
}
},
"if1": {
"process_id": "if",
"arguments": {
"accept": 1,
"reject": {
"from_node": "power1"
},
"value": {
"from_node": "isnodata1"
}
}
},
"arrayelement3": {
"process_id": "array_element",
"arguments": {
"data": {
"from_parameter": "data"
},
"index": 1
}
},
"log2": {
"process_id": "log",
"arguments": {
"base": 10,
"x": {
"from_node": "arrayelement3"
}
}
},
"multiply2": {
"process_id": "multiply",
"arguments": {
"x": 10,
"y": {
"from_node": "log2"
}
}
},
"add2": {
"process_id": "add",
"arguments": {
"x": {
"from_node": "multiply2"
},
"y": 83
}
},
"divide2": {
"process_id": "divide",
"arguments": {
"x": {
"from_node": "add2"
},
"y": 20
}
},
"power2": {
"process_id": "power",
"arguments": {
"base": 10,
"p": {
"from_node": "divide2"
}
}
},
"arrayelement4": {
"process_id": "array_element",
"arguments": {
"data": {
"from_parameter": "data"
},
"index": 1
}
},
"isnodata2": {
"process_id": "is_nodata",
"arguments": {
"x": {
"from_node": "arrayelement4"
}
}
},
"if2": {
"process_id": "if",
"arguments": {
"accept": 1,
"reject": {
"from_node": "power2"
},
"value": {
"from_node": "isnodata2"
}
}
},
"arraycreate1": {
"process_id": "array_create",
"arguments": {
"data": [
{
"from_node": "if1"
},
{
"from_node": "if2"
}
]
},
"result": true
}
}
}
}
},
"apply2": {
"process_id": "apply",
"arguments": {
"data": {
"from_node": "applydimension3"
},
"process": {
"process_graph": {
"linearscalerange2": {
"process_id": "linear_scale_range",
"arguments": {
"inputMax": 65534,
"inputMin": 1,
"outputMax": 65534,
"outputMin": 1,
"x": {
"from_parameter": "x"
}
},
"result": true
}
}
}
}
},
"resamplecubespatial1": {
"process_id": "resample_cube_spatial",
"arguments": {
"data": {
"from_node": "apply1"
},
"method": "near",
"target": {
"from_node": "apply2"
}
}
},
"mergecubes1": {
"process_id": "merge_cubes",
"arguments": {
"cube1": {
"from_node": "resamplecubespatial1"
},
"cube2": {
"from_node": "apply2"
}
}
},
"filtertemporal1": {
"process_id": "filter_temporal",
"arguments": {
"data": {
"from_node": "mergecubes1"
},
"extent": [
"2018-03-01",
"2018-10-31"
]
}
},
"loadcollection3": {
"process_id": "load_collection",
"arguments": {
"id": "COPERNICUS_30",
"spatial_extent": {
"east": 4660000,
"south": 3240000,
"west": 4640000,
"north": 3260000,
"crs": "EPSG:3035"
},
"temporal_extent": null
}
},
"resamplespatial3": {
"process_id": "resample_spatial",
"arguments": {
"align": "upper-left",
"data": {
"from_node": "loadcollection3"
},
"method": "cubic",
"projection": 3035,
"resolution": 10
}
},
"reducedimension1": {
"process_id": "reduce_dimension",
"arguments": {
"data": {
"from_node": "resamplespatial3"
},
"dimension": "t",
"reducer": {
"process_graph": {
"max1": {
"process_id": "max",
"arguments": {
"data": {
"from_parameter": "data"
}
},
"result": true
}
}
}
}
},
"mergecubes2": {
"process_id": "merge_cubes",
"arguments": {
"cube1": {
"from_node": "filtertemporal1"
},
"cube2": {
"from_node": "reducedimension1"
}
}
},
"apply3": {
"process_id": "apply",
"arguments": {
"data": {
"from_node": "mergecubes2"
},
"process": {
"process_graph": {
"linearscalerange3": {
"process_id": "linear_scale_range",
"arguments": {
"inputMax": 65534,
"inputMin": 0,
"outputMax": 65534,
"outputMin": 0,
"x": {
"from_parameter": "x"
}
},
"result": true
}
}
}
}
},
"applyneighborhood1": {
"process_id": "apply_neighborhood",
"arguments": {
"data": {
"from_node": "apply3"
},
"overlap": [],
"process": {
"process_graph": {
"runudf1": {
"process_id": "run_udf",
"arguments": {
"context": {
"start_month": 3,
"end_month": 10,
"year": 2018,
"modeldir": "tmp/model/",
"modeltag": "20230615T144208-24ts-hrlvlcc-v200",
"target_crs": 3035,
"s1_orbitdirection": "DESCENDING",
"apply_bvl_mask": false,
"all_probabilities": true,
"provider": "creodias",
"tiles": "33UV*",
"METEO_data": {
"dims": [
"t",
"bands",
"x",
"y"
],
"attrs": {
"nodata": 65535,
"crs": "+proj=longlat +datum=WGS84 +no_defs ",
"grid_mapping": "spatial_ref",
"dtype": "uint16",
"shape": [
36,
1,
4,
2
]
},
"data": [
[
[
[
27657,
27658
],
[
27661,
27656
],
[
27654,
27638
],
[
27644,
27621
]
]
],
[
[
[
27372,
27374
],
[
27365,
27360
],
[
27356,
27339
],
[
27344,
27321
]
]
],
[
[
[
27760,
27762
],
[
27740,
27736
],
[
27736,
27719
],
[
27726,
27703
]
]
],
[
[
[
27194,
27197
],
[
27206,
27201
],
[
27196,
27179
],
[
27184,
27160
]
]
],
[
[
[
27321,
27323
],
[
27327,
27321
],
[
27318,
27300
],
[
27307,
27282
]
]
],
[
[
[
26740,
26744
],
[
26734,
26727
],
[
26720,
26701
],
[
26706,
26678
]
]
],
[
[
[
27172,
27175
],
[
27161,
27154
],
[
27151,
27131
],
[
27139,
27112
]
]
],
[
[
[
27520,
27522
],
[
27534,
27528
],
[
27528,
27507
],
[
27518,
27489
]
]
],
[
[
[
27610,
27612
],
[
27616,
27610
],
[
27611,
27590
],
[
27601,
27572
]
]
],
[
[
[
28354,
28355
],
[
28361,
28356
],
[
28363,
28341
],
[
28357,
28327
]
]
],
[
[
[
28774,
28774
],
[
28780,
28774
],
[
28785,
28763
],
[
28782,
28750
]
]
],
[
[
[
28720,
28720
],
[
28732,
28726
],
[
28736,
28714
],
[
28732,
28700
]
]
],
[
[
[
28809,
28809
],
[
28824,
28818
],
[
28828,
28805
],
[
28825,
28792
]
]
],
[
[
[
28912,
28912
],
[
28929,
28922
],
[
28934,
28910
],
[
28930,
28896
]
]
],
[
[
[
29424,
29422
],
[
29433,
29427
],
[
29442,
29418
],
[
29441,
29406
]
]
],
[
[
[
29416,
29414
],
[
29424,
29416
],
[
29431,
29406
],
[
29430,
29395
]
]
],
[
[
[
29255,
29253
],
[
29264,
29255
],
[
29268,
29243
],
[
29265,
29230
]
]
],
[
[
[
29050,
29048
],
[
29055,
29046
],
[
29057,
29031
],
[
29052,
29017
]
]
],
[
[
[
29251,
29249
],
[
29258,
29249
],
[
29260,
29235
],
[
29257,
29222
]
]
],
[
[
[
29262,
29259
],
[
29263,
29254
],
[
29264,
29239
],
[
29261,
29226
]
]
],
[
[
[
29700,
29696
],
[
29690,
29681
],
[
29695,
29671
],
[
29693,
29660
]
]
],
[
[
[
29842,
29839
],
[
29848,
29840
],
[
29855,
29831
],
[
29854,
29821
]
]
],
[
[
[
29470,
29467
],
[
29475,
29466
],
[
29478,
29454
],
[
29474,
29442
]
]
],
[
[
[
29177,
29174
],
[
29185,
29175
],
[
29184,
29160
],
[
29178,
29147
]
]
],
[
[
[
29175,
29172
],
[
29184,
29175
],
[
29183,
29161
],
[
29178,
29148
]
]
],
[
[
[
29167,
29164
],
[
29185,
29177
],
[
29184,
29163
],
[
29179,
29150
]
]
],
[
[
[
28523,
28520
],
[
28527,
28517
],
[
28519,
28498
],
[
28511,
28482
]
]
],
[
[
[
28474,
28472
],
[
28485,
28475
],
[
28477,
28457
],
[
28468,
28441
]
]
],
[
[
[
28731,
28729
],
[
28744,
28736
],
[
28740,
28721
],
[
28733,
28707
]
]
],
[
[
[
28171,
28170
],
[
28171,
28163
],
[
28162,
28143
],
[
28152,
28127
]
]
],
[
[
[
28272,
28271
],
[
28290,
28282
],
[
28283,
28265
],
[
28274,
28250
]
]
],
[
[
[
27841,
27840
],
[
27839,
27832
],
[
27829,
27811
],
[
27817,
27794
]
]
],
[
[
[
27431,
27432
],
[
27423,
27415
],
[
27409,
27391
],
[
27396,
27372
]
]
],
[
[
[
27862,
27862
],
[
27858,
27851
],
[
27849,
27833
],
[
27839,
27816
]
]
],
[
[
[
27431,
27432
],
[
27424,
27418
],
[
27413,
27396
],
[
27400,
27378
]
]
],
[
[
[
27804,
27805
],
[
27794,
27789
],
[
27787,
27771
],
[
27777,
27755
]
]
]
],
"coords": {
"bands": {
"dims": [
"bands"
],
"attrs": {
"dtype": "<U16",
"shape": [
1
]
},
"data": [
"temperature_mean"
]
},
"x": {
"dims": [
"x"
],
"attrs": {
"dtype": "float64",
"shape": [
4
]
},
"data": [
14.717024556207523,
14.817024556207523,
14.917024556207522,
15.017024556207522
]
},
"y": {
"dims": [
"y"
],
"attrs": {
"dtype": "float64",
"shape": [
2
]
},
"data": [
52.207046686641945,
52.307046686641954
]
},
"t": {
"dims": [
"t"
],
"attrs": {
"dtype": "datetime64[ns]",
"shape": [
36
]
},
"data": [
"2018-01-01 00:00:00",
"2018-01-11 00:00:00",
"2018-01-21 00:00:00",
"2018-02-01 00:00:00",
"2018-02-11 00:00:00",
"2018-02-21 00:00:00",
"2018-03-01 00:00:00",
"2018-03-11 00:00:00",
"2018-03-21 00:00:00",
"2018-04-01 00:00:00",
"2018-04-11 00:00:00",
"2018-04-21 00:00:00",
"2018-05-01 00:00:00",
"2018-05-11 00:00:00",
"2018-05-21 00:00:00",
"2018-06-01 00:00:00",
"2018-06-11 00:00:00",
"2018-06-21 00:00:00",
"2018-07-01 00:00:00",
"2018-07-11 00:00:00",
"2018-07-21 00:00:00",
"2018-08-01 00:00:00",
"2018-08-11 00:00:00",
"2018-08-21 00:00:00",
"2018-09-01 00:00:00",
"2018-09-11 00:00:00",
"2018-09-21 00:00:00",
"2018-10-01 00:00:00",
"2018-10-11 00:00:00",
"2018-10-21 00:00:00",
"2018-11-01 00:00:00",
"2018-11-11 00:00:00",
"2018-11-21 00:00:00",
"2018-12-01 00:00:00",
"2018-12-11 00:00:00",
"2018-12-21 00:00:00"
]
},
"spatial_ref": {
"dims": [],
"attrs": {
"spatial_ref": "GEOGCRS[\"unknown\",DATUM[\"World Geodetic System 1984\",ELLIPSOID[\"WGS 84\",6378137,298.257223563,LENGTHUNIT[\"metre\",1]],ID[\"EPSG\",6326]],PRIMEM[\"Greenwich\",0,ANGLEUNIT[\"degree\",0.0174532925199433],ID[\"EPSG\",8901]],CS[ellipsoidal,2],AXIS[\"longitude\",east,ORDER[1],ANGLEUNIT[\"degree\",0.0174532925199433,ID[\"EPSG\",9122]]],AXIS[\"latitude\",north,ORDER[2],ANGLEUNIT[\"degree\",0.0174532925199433,ID[\"EPSG\",9122]]]]",
"crs_wkt": "GEOGCRS[\"unknown\",DATUM[\"World Geodetic System 1984\",ELLIPSOID[\"WGS 84\",6378137,298.257223563,LENGTHUNIT[\"metre\",1]],ID[\"EPSG\",6326]],PRIMEM[\"Greenwich\",0,ANGLEUNIT[\"degree\",0.0174532925199433],ID[\"EPSG\",8901]],CS[ellipsoidal,2],AXIS[\"longitude\",east,ORDER[1],ANGLEUNIT[\"degree\",0.0174532925199433,ID[\"EPSG\",9122]]],AXIS[\"latitude\",north,ORDER[2],ANGLEUNIT[\"degree\",0.0174532925199433,ID[\"EPSG\",9122]]]]",
"dtype": "int64",
"shape": []
},
"data": 0
}
},
"name": null
}
},
"data": {
"from_parameter": "data"
},
"runtime": "Python",
"udf": "# -*- coding: utf-8 -*-\n# Uncomment the import only for coding support\nfrom openeo.udf import XarrayDataCube\nfrom typing import Dict\nimport xarray as xr\nimport numpy as np\nimport threading\nimport functools\nimport datetime\nfrom pathlib import Path\nimport sys\nimport pandas as pd\n\nfrom typing import List\n\n_threadlocal = threading.local()\n\nDEFAULTMODELDIR = Path('/vitodata/EEA_HRL_VLCC/data/inference/crop_type/models/') # NOQA\nDEFAULTMODELTAG = '20221102T184737-transformer_optical_dem'\n\n\ndef _setup_logging():\n global logger\n try:\n from loguru import logger\n except ImportError:\n import logging\n logger = logging.getLogger(__name__)\n\n\ndef scaled_s1_to_pwr(scaled_db):\n from cropclass.features import to_pwr\n '''Method to remove the scaling from S1 data\n fetched from openeo preprocessing workflow\n and converting to power\n '''\n\n result = to_pwr(20 * np.log10(scaled_db) - 83)\n result[~np.isfinite(result)] = np.nan\n\n return result\n\n\ndef get_preproc_settings():\n from cropclass.features import divide_100, divide_10000, to_db # NOQA\n return {\n 'S2':\n {\n 'dtype': np.float32,\n 'src_nodata': 65535, # Nodata value for UINT16 in OpenEO\n 'dst_nodata': np.nan,\n 'pre_func': divide_10000,\n 'interpolate': True\n },\n 'S1':\n {\n 'dtype': np.float32,\n 'src_nodata': 65535, # Nodata value for UINT16 in OpenEO\n 'dst_nodata': np.nan,\n 'bands': ['VV', 'VH'],\n 'pre_func': scaled_s1_to_pwr,\n 'post_func': to_db,\n 'interpolate': True\n },\n 'METEO':\n {\n 'dtype': np.float32,\n 'src_nodata': 65535, # Nodata value for UINT16 in OpenEO\n 'dst_nodata': np.nan,\n 'pre_func': divide_100,\n 'interpolate': False # Should not be needed\n },\n }\n\n\ndef convert_meteojson_to_xarray(xrarray_frame, meteo_json, startdate, enddate):\n \"\"\"\n gets all data from meteo_json temoral filter this data to startdata\n and endate. Calculates the spatial average (with nodata values masked.)\n and pushes it into a xr-array.\n returns the data in the xarray\n @param xrarray_frame: xarray frame which should be have the correct coords\n @param meteo_json: input json with the meteo data from agera5 as extracted\n by preprocessing_agera5 script\n @param startdate: the startdate of the data that is needed\n @param enddate: the enddate of the data that is needed\n @return: xarray\n \"\"\"\n\n # aggregate data\n meteo_data = meteo_json.get(\"data\")\n dates = [datetime.datetime.strptime(date, \"%Y-%m-%d %H:%M:%S\")\n for date in meteo_json.get(\"coords\").get(\"t\").get(\"data\")]\n masked_meteo = np.ma.masked_values(\n meteo_data, meteo_json.get(\"attrs\").get(\"nodata\"))\n aggregate_meteo = masked_meteo.mean(axis=(1, 2, 3), keepdims=False)\n aggregate_meteo.filled(fill_value=meteo_json.get(\"attrs\").get(\"nodata\"))\n # filter based on timings\n enddate = pd.to_datetime(enddate)\n startdate = pd.to_datetime(startdate)\n filter = [date >= startdate and date < enddate for date in dates]\n aggregate_meteo = aggregate_meteo[filter]\n # tile and assign to xarray\n res_agg = aggregate_meteo.reshape((1, len(aggregate_meteo), 1, 1))\n tiled_data = np.tile(\n res_agg, (1, 1, len(xrarray_frame['x']), len(xrarray_frame['y'])))\n xarray = xr.ones_like(xrarray_frame)\n xarray.data = tiled_data\n xarray = xarray.assign_coords(\n {'bands': meteo_json.get(\"coords\").get(\"bands\").get(\"data\")})\n\n return xarray\n\n\ndef add_dem_features(features, dem_data):\n from satio.utils import dem_attrs\n from satio.features import Features\n from cropclass.features.inference import UnexpectedRangeError\n\n dem_data = dem_data.astype(np.float32)\n dem_data[dem_data < -10000] = np.nan\n\n MAX_ALLOWED_ALTITUDE = 8000\n MIN_ALLOWED_ALTITUDE = -50\n\n # Some checks on valid DEM data\n if np.nanmax(dem_data) > MAX_ALLOWED_ALTITUDE:\n raise UnexpectedRangeError(\n ('Detected unrealistic DEM data: '\n f'{np.nanmax(dem_data)} is higher '\n f'than max allowed {MAX_ALLOWED_ALTITUDE}'))\n if np.nanmin(dem_data) < MIN_ALLOWED_ALTITUDE:\n raise UnexpectedRangeError(\n ('Detected unrealistic DEM data: '\n f'{np.nanmin(dem_data)} is lower '\n f'than min allowed {MIN_ALLOWED_ALTITUDE}'))\n\n # To compute slope, nodata value is expected to be -9999\n slope, _ = dem_attrs(np.where(np.isfinite(dem_data), dem_data, -9999))\n slope[slope == -9999] = np.nan\n\n dem_arr = np.array([dem_data, slope])\n dem_feats = Features(dem_arr, ['DEM-alt-20m', 'DEM-slo-20m'])\n\n return features.merge(dem_feats)\n\n\n@functools.lru_cache(maxsize=25)\ndef load_croptype_model(modeldir: str, modeltag: str):\n\n from vito_crop_classification.model import Model\n from cropclass.utils import download_and_unpack\n\n if str(modeldir).startswith('http'):\n modeldir, modeltag = download_and_unpack(modeldir, modeltag)\n\n # load a model from its tag\n model = Model.load(mdl_f=Path(modeldir) / modeltag)\n\n # Get classes that will be mapped\n class_names = model.get_class_names()\n\n return model, class_names\n\n\ndef _make_prediction(features, modeldir, modeltag,\n apply_worldcovermask=True, ignore_classes=[]):\n from vito_crop_classification import run_prediction\n\n model, classes = load_croptype_model(modeldir, modeltag)\n\n df = features.df\n valid_indices = df[~df.isna().any(axis=1)].index\n #valid_indices = pd.Index([2,4,5,6,7,10,100,101,102])\n #invalid_indices = df[df.isna().any(axis=1)].index\n df_filter = df.loc[df.index.isin(valid_indices), slice(None)]\n df_filter.reset_index(inplace = True)\n\n\n pred_columns = ['prediction_id', 'prediction_name', 'probability',\n 'probabilities', 'embedding']\n\n df = df.reindex(columns=[*df.columns,*pred_columns])\n\n # initialize predictor over new data\n logger.info('Running Predictor.')\n import torch\n with torch.no_grad():\n if len(valid_indices) != 0:\n df.loc[df.index.isin(valid_indices), pred_columns] = run_prediction(\n df=df_filter,\n model=model,\n patch_smoothing=False,\n transform=True,\n ignore_classes=ignore_classes).set_index(valid_indices)\n\n # Set NaN strings where we don't have prediction for later\n # translation to raster values\n df.loc[~df.index.isin(valid_indices),\n ['prediction_id', 'prediction_name']] = 'NaN'\n df.loc[~df.index.isin(valid_indices), 'probability'] = 0\n df.loc[~df.index.isin(valid_indices),\n 'probabilities'] = (\n df.loc[~df.index.isin(valid_indices), 'probabilities'].apply(\n lambda _: np.zeros(len(classes), dtype=np.float32)))\n\n # Get prediction and probability arrays\n prediction = df['prediction_id'].values\n probability = df['probability'].values\n probabilities = df['probabilities'].values\n\n # Rounding to the closest integer, obtaining a value between\n # 0 and 100 encoded in bytes\n probability = ((probability * 100.0) + 0.5).astype(np.uint8)\n probabilities = np.moveaxis(\n (np.stack(probabilities) * 100.0 + 0.5).astype(np.uint8), 0, 1\n )\n\n if apply_worldcovermask:\n logger.info('Applying WorldCover mask.')\n worldcoverlabel = features.df['WORLDCOVER-label-10m']\n prediction[~(worldcoverlabel.isin([30, 40]))] = \"NoCrop\"\n\n return (\n prediction.reshape(1, *features.data.shape[1:3]),\n probability.reshape(1, *features.data.shape[1:3]),\n probabilities.reshape(len(classes), *features.data.shape[1:3]),\n classes\n )\n\n\ndef _compute_features(inarr: xr.DataArray,\n startdate: str,\n enddate: str,\n n_tsteps: int,\n apply_worldcovermask: bool = False,\n segmentation: bool = False,\n segm_settings: Dict = None,\n outofrange_toleranceratio: np.float32 = 0.001,\n **kwargs) -> xr.DataArray:\n from cropclass.features import FeaturesComputer\n from cropclass.features.settings import get_feature_settings\n from satio.features import Features\n\n # Fix to see if the logger is defined or not, as the compute_features\n # function can be called alone\n try:\n logger\n except NameError:\n _setup_logging()\n\n # Make sure we have a controlled dimension order\n inarr = inarr.transpose('bands', 't', 'x', 'y')\n\n # Add meteo band if needed\n bands = list(inarr.bands.values)\n if 'METEO_data' in kwargs:\n if kwargs.get('METEO_data') is not None:\n meteo_array = convert_meteojson_to_xarray(\n inarr.sel(bands=[bands[0]]), kwargs.get('METEO_data'),\n startdate, enddate)\n inarr = xr.concat([inarr, meteo_array], dim=\"bands\")\n\n # Initialize features computer\n fc = FeaturesComputer(\n get_feature_settings(n_tsteps),\n startdate,\n enddate,\n preprocessing_settings=get_preproc_settings(),\n segmentation=segmentation,\n segm_settings=segm_settings,\n outofrange_toleranceratio=outofrange_toleranceratio\n )\n\n # Identify the S1, S2 and METEO bands\n bands = list(inarr.bands.values)\n S2_bands = [b for b in bands if 'B' in b]\n S1_bands = [b for b in bands if 'V' in b]\n\n # Only the temperature-mean band is useful now\n METEO_bands = [b for b in bands if 'temperature_mean' in b]\n\n # Construct sensor-specific arrays.\n sensor_arrays = {'S2': inarr.sel(bands=S2_bands)}\n\n # Add S1 features, if exists\n # at this point they will still be scaled!\n if (len(S1_bands)) > 0:\n sensor_arrays['S1'] = inarr.sel(bands=S1_bands)\n\n # Add Meteo features, if exists\n if len(METEO_bands) > 0:\n sensor_arrays['METEO'] = inarr.sel(bands=METEO_bands)\n\n # Compute features\n features = fc.get_features(sensor_arrays)\n\n # Add DEM\n features = add_dem_features(features, inarr.sel(\n bands='DEM').max(dim='t').values)\n\n # Add worldcover\n if apply_worldcovermask:\n worldcover_feats = Features(inarr.sel(bands='MAP').max(\n dim='t').values, ['WORLDCOVER-label-10m'])\n features = features.merge(worldcover_feats)\n\n return features\n\n\ndef classify(\n inarr: xr.DataArray,\n start_month: int,\n end_month: int,\n year: int,\n modeltag: str = DEFAULTMODELTAG,\n modeldir: Path = DEFAULTMODELDIR,\n apply_worldcovermask: bool = False,\n ignore_classes: List[str] = [],\n all_probabilities: bool = False,\n **processing_opts) -> xr.DataArray:\n \"\"\"Main inference classification UDF for OpenEO workflow.\n This function takes an input xarray fetched through OpenEO\n and returns a result xarray back to OpenEO for saving.\n\n Args:\n inarr (xr.DataArray): preprocessed\n input arrays for the different sensors\n start_month (int): start month for the feature computation\n end_month (int): end month for the feature computation\n year (int): Year in which mapped season ends.\n apply_worldcovermask (bool, optional): Whether or not to apply\n mask based on worldcover 30 and 40 values (grass and crop).\n Pixels outside this mask will get value 255 in\n output product. Defaults to False.\n ignore_classes (str, optional): pass list of classes to ignore\n custom_dependency_path (str, optional: optional path to be\n added in front of the python path\n\n Returns:\n xr.DataArray: resulting classification xarray.\n \"\"\"\n import cropclass\n from cropclass.utils.seasons import get_processing_dates\n from cropclass.features.settings import compute_n_tsteps\n import vito_crop_classification\n\n try:\n logger\n except NameError:\n _setup_logging()\n\n logger.info(f'Using cropclass library from: {cropclass.__file__}')\n logger.info(('Using vito_crop_classification '\n f'library from: {vito_crop_classification.__file__}'))\n\n # Infer exact processing dates from months and year\n startdate, enddate = get_processing_dates(start_month, end_month, year)\n logger.info(f'Inferred processing date: {startdate} - {enddate}')\n\n # Infer required amount of tsteps\n n_tsteps = compute_n_tsteps(start_month, end_month)\n logger.info(f'Inferred n_tsteps: {n_tsteps}')\n\n # Store the original dimension order for later\n orig_dims = list(inarr.dims)\n orig_dims.remove('t') # Time dimension will be removed after the inference\n\n features = _compute_features(inarr, startdate, enddate, n_tsteps,\n apply_worldcovermask, **processing_opts)\n\n # Make the crop type prediction\n logger.info(f'Using model: {modeltag}')\n prediction, probability, cls_probabilities, cls_names = _make_prediction(\n features,\n modeldir,\n modeltag=modeltag,\n apply_worldcovermask=apply_worldcovermask,\n ignore_classes=ignore_classes\n )\n\n features = None\n\n if not all_probabilities:\n # Zip the prediction and probability arrays together\n predicted_data = np.array([prediction.squeeze(),\n probability.squeeze()])\n\n # Finally transform result to DataArray\n result_da = xr.DataArray(\n predicted_data,\n coords=[\n ['croptype', 'probability'],\n inarr.coords[\"x\"],\n inarr.coords[\"y\"]\n ],\n dims=[\"bands\", \"x\", \"y\"]\n )\n else:\n # Adds the prediciton and winning proba bands, then the class\n # probabilities one by one\n band_list = [prediction.squeeze(), probability.squeeze()]\n for prob_idx in range(cls_probabilities.shape[0]):\n band_list.append(cls_probabilities[prob_idx].squeeze())\n band_coords = ['croptype', 'probability']\n for cls_name in cls_names:\n band_coords.append('probability_' +\n cls_name.replace(' ', '_').replace(',', ''))\n\n predicted_data = np.array(band_list)\n\n result_da = xr.DataArray(\n predicted_data,\n coords=[\n band_coords,\n inarr.coords['x'],\n inarr.coords['y']\n ],\n dims=['bands', 'x', 'y']\n )\n\n # result_da = stacked_array.assign_coords(\n # {\"feature\": ft_labels_final}).rename(feature='bands')\n\n # And make sure we revert back to original dimension order\n result_da = result_da.transpose(*orig_dims)\n\n return result_da\n\n\ndef apply_datacube(cube: XarrayDataCube, context: Dict) -> XarrayDataCube:\n sys.path.insert(0, 'tmp/venv_static')\n sys.path.insert(0, 'tmp/cropclasslib')\n sys.path.insert(0, 'tmp/vitocropclassification')\n\n from cropclass.postprocessing.layer_format import translate_layer\n\n _setup_logging()\n\n # Extract xarray.DataArray from the cube\n inarr = cube.get_array()\n\n # Run the two-stage classification\n predictions = classify(inarr, **context)\n\n translated_predictions = translate_layer(\n predictions, translation_table=context.get('translation_table', None))\n\n # Wrap result in an OpenEO datacube\n return XarrayDataCube(translated_predictions)\n"
},
"result": true
}
}
},
"size": [
{
"dimension": "x",
"value": 128,
"unit": "px"
},
{
"dimension": "y",
"value": 128,
"unit": "px"
}
]
}
},
"apply4": {
"process_id": "apply",
"arguments": {
"data": {
"from_node": "applyneighborhood1"
},
"process": {
"process_graph": {
"linearscalerange4": {
"process_id": "linear_scale_range",
"arguments": {
"inputMax": 10000,
"inputMin": 0,
"outputMax": 10000,
"outputMin": 0,
"x": {
"from_parameter": "x"
}
},
"result": true
}
}
}
}
},
"dropdimension1": {
"process_id": "drop_dimension",
"arguments": {
"data": {
"from_node": "apply4"
},
"name": "t"
}
},
"renamelabels1": {
"process_id": "rename_labels",
"arguments": {
"data": {
"from_node": "dropdimension1"
},
"dimension": "bands",
"target": [
"croptype",
"probability"
]
}
},
"saveresult1": {
"process_id": "save_result",
"arguments": {
"data": {
"from_node": "renamelabels1"
},
"format": "GTiff",
"options": {
"overview_method": "mode",
"filename_prefix": "CROP_2018_E464N324-03035-010m_V100"
}
},
"result": true
}
}
}
When ignoring a tile, it could be replaced with NODATA. However, the
SourceNamecan beStringNamethen, and would cause issues in the FIleLayerProviderj-231219c859b44e60a252654db30210bbSuppressed: scala.MatchError: (GridBoundsRasterRegion(MosaicRasterSource(List(org.openeo.geotrellis.layers.NoDataRasterSource(uint8, GridExtent(Extent(4561530.0, 3144110.0, 4678950.0, 3261370.0), CellSize(10.0,10.0), 11742x11726), EPSG:3035), org.openeo.geotrellis.layers.NoDataRasterSource(uint8, GridExtent(Extent(4561530.0, 3144110.0, 4678950.0, 3261370.0), CellSize(10.0,10.0), 11742x11726), EPSG:3035), org.openeo.geotrellis.layers.NoDataRasterSource(uint8, GridExtent(Extent(4561530.0, 3144110.0, 4678950.0, 3261370.0), CellSize(10.0,10.0), 11742x11726), EPSG:3035), org.openeo.geotrellis.layers.ValueOffsetRasterSource@5e9db88c, org.openeo.geotrellis.layers.ValueOffsetRasterSource@1899243b, org.openeo.geotrellis.layers.ValueOffsetRasterSource@62fb1f01, org.openeo.geotrellis.layers.NoDataRasterSource(uint8, GridExtent(Extent(4561530.0, 3144110.0, 4678950.0, 3261370.0), CellSize(10.0,10.0), 11742x11726), EPSG:3035), org.openeo.geotrellis.layers.ValueOffsetRasterSource@733bbf7f, org.openeo.geotrellis.layers.ValueOffsetRasterSource@65c60562, GDALRasterSource(/vsis3/EODATA/Sentinel-2/MSI/L2A/2018/10/28/S2B_MSIL2A_20181028T100119_N0206_R122_T33UVT_20181028T135751.SAFE/GRANULE/L2A_T33UVT_A008585_20181028T100117/IMG_DATA/R20m/T33UVT_20181028T100119_SCL_20m.jp2,GDALWarpOptions(-of VRT -r near -tap -tr 10.0 10.0 -t_srs +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs -ovr AUTO -te 4561530.0 3144110.0 4678950.0 3261370.0 -te_srs +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs ))), EPSG:3035, GridExtent(Extent(4561530.0, 3144110.0, 4678950.0, 3261370.0), CellSize(10.0,10.0), 11742x11726), StringName(org.openeo.geotrellis.layers.NoDataRasterSource(uint8, GridExtent(Extent(4561530.0, 3144110.0, 4678950.0, 3261370.0), CellSize(10.0,10.0), 11742x11726), EPSG:3035))),GridBounds(9383,137,9894,648)),StringName(org.openeo.geotrellis.layers.NoDataRasterSource(uint8, GridExtent(Extent(4561530.0, 3144110.0, 4678950.0, 3261370.0), CellSize(10.0,10.0), 11742x11726), EPSG:3035))) (of class scala.Tuple2)