# import libraries
import os
import sys
import ilwis
import numpy as np
from netCDF4 import Dataset
from osgeo import gdal, osr
from osgeo.gdalconst import *
import matplotlib.pyplot as plt
import openeo
openeo.client_version()
'0.25.0'
!rm -rf /home/eoafrica/.local/share/openeo-python-client/refresh-tokens.json
con = openeo.connect(url="openeo.vito.be")
con.authenticate_oidc()
Authenticated using refresh token.
<Connection to 'https://openeo.vito.be/openeo/1.2/' with OidcBearerAuth>
print(con.list_collection_ids())
['MAPEO_WATER_TUR_V1', 'COP_DEM_EU_25M', 'ESA_WORLDCEREAL_ACTIVECROPLAND', 'ESA_WORLDCEREAL_IRRIGATION', 'ESA_WORLDCEREAL_TEMPORARYCROPS', 'ESA_WORLDCEREAL_WINTERCEREALS', 'ESA_WORLDCEREAL_MAIZE', 'ESA_WORLDCEREAL_SPRINGCEREALS', 'SENTINEL1_GRD_SIGMA0', 'S1_GRD_SIGMA0_ASCENDING', 'S1_GRD_SIGMA0_DESCENDING', 'SENTINEL3_SYNERGY_VG1', 'SENTINEL3_SYNERGY_VG10', 'TERRASCOPE_S2_FAPAR_V2', 'TERRASCOPE_S2_NDVI_V2', 'TERRASCOPE_S2_LAI_V2', 'TERRASCOPE_S2_FCOVER_V2', 'TERRASCOPE_S2_TOC_V2', 'TERRASCOPE_S1_SLC_COHERENCE_V1', 'SENTINEL1_GAMMA0_SENTINELHUB', 'SENTINEL1_GRD', 'SENTINEL2_L1C_SENTINELHUB', 'SENTINEL2_L2A_SENTINELHUB', 'SENTINEL2_L2A_MOSAIC_120', 'PROBAV_L3_S10_TOC_333M', 'PROBAV_L3_S5_TOC_100M', 'PROBAV_L3_S1_TOC_100M', 'PROBAV_L3_S1_TOC_333M', 'TERRASCOPE_S5P_L3_CO_TD_V1', 'TERRASCOPE_S5P_L3_CO_TM_V1', 'TERRASCOPE_S5P_L3_CO_TY_V1', 'LANDSAT8-9_L1', 'LANDSAT8-9_L2', 'LANDSAT4-5_TM_L1', 'LANDSAT4-5_TM_L2', 'LANDSAT7_ETM_L1', 'LANDSAT7_ETM_L2', 'MODIS', 'SENTINEL3_OLCI_L1B', 'SENTINEL3_SLSTR', 'COPERNICUS_30', 'COPERNICUS_90', 'TERRASCOPE_S2_RHOW_V1', 'TERRASCOPE_S2_CHL_V1', 'TERRASCOPE_S2_SPM_V1', 'TERRASCOPE_S2_TUR_V1', 'ESA_WORLDCOVER_10M_2020_V1', 'ESA_WORLDCOVER_10M_2021_V2', 'CORINE_LAND_COVER', 'CORINE_LAND_COVER_ACCOUNTING_LAYERS', 'GHS_BUILT_S2', 'GLOBAL_LAND_COVER', 'GLOBAL_SURFACE_WATER', 'SEA_ICE_INDEX', 'SENTINEL_5P_L2', 'MAPZEN_DEM', 'LANDSAT1-5_MSS_L1', 'WATER_BODIES', 'SENTINEL1_CARD4L', 'VEGETATION_INDICES', 'SEASONAL_TRAJECTORIES', 'VEGETATION_PHENOLOGY_AND_PRODUCTIVITY_PARAMETERS_SEASON_1', 'VEGETATION_PHENOLOGY_AND_PRODUCTIVITY_PARAMETERS_SEASON_2', 'FRACTIONAL_SNOW_COVER', 'OPENEO_CROPTYPE_2021', 'EGM2008', 'CGLS_NDVI_LTS_V2_GLOBAL', 'CGLS_GDMP300_V1_GLOBAL', 'CGLS_GDMP_V2_GLOBAL', 'CGLS_SSM_V1_EUROPE', 'CGLS_FAPAR_V2_GLOBAL', 'CGLS_LAI_V2_GLOBAL', 'CGLS_FCOVER_V2_GLOBAL', 'CGLS_LAI300_V1_GLOBAL', 'CGLS_NDVI300_V1_GLOBAL', 'CGLS_FCOVER300_V1_GLOBAL', 'CGLS_FAPAR300_V1_GLOBAL', 'CGLS_NDVI300_V2_GLOBAL', 'CGLS_NDVI_V3_GLOBAL', 'CGLS_NDVI_V2_GLOBAL', 'CGLS_SWI_V1_EUROPE', 'TERRASCOPE_S1_GAMMA0_V1', 'SENTINEL2_L1C_INCD', 'SENTINEL2_L2A', 'TERRASCOPE_S5P_L3_NO2_TD_V1', 'TERRASCOPE_S5P_L3_NO2_TM_V1', 'TERRASCOPE_S5P_L3_NO2_TY_V1', 'AGERA5', 'PLANETSCOPE', 'POPULATION_DENSITY', 'SENTINEL2_L1C']
con.describe_collection('GLOBAL_LAND_COVER')#web listing
collection_id = 'GLOBAL_LAND_COVER'
spatial_extent = {"west": 32.37, "south": -20.941, "east": 34.758, "north": -18.669}
cube = con.load_collection(
"GLOBAL_LAND_COVER",
spatial_extent=spatial_extent,
temporal_extent=["2015-01-01", "2019-12-31"],
bands=["Discrete_Classification_map"],
)
%%time
cube.download('data/Buzi_LCC.nc')
input_file = (r'data/Buzi_LCC.nc')
print(input_file)
data/Buzi_LCC.nc
dataset = gdal.Open(input_file)
nc_info = gdal.Info(dataset)
print (nc_info)
Driver: netCDF/Network Common Data Format Files: data/Buzi_LCC.nc Size is 2580, 2454 Coordinate System is: GEOGCRS["WGS 84", 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]], AXIS["latitude",north, ORDER[2], ANGLEUNIT["degree",0.0174532925199433]]] Data axis to CRS axis mapping: 1,2 Origin = (32.369999999999997,-18.669000000000000) Pixel Size = (0.000925925925926,-0.000925925925926) Metadata: crs#crs_wkt=GEOGCS["WGS 84", DATUM["World Geodetic System 1984", SPHEROID["WGS 84", 6378137.0, 298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich", 0.0, AUTHORITY["EPSG","8901"]], UNIT["degree", 0.017453292519943295], AXIS["Geodetic longitude", EAST], AXIS["Geodetic latitude", NORTH], AUTHORITY["EPSG","4326"]] crs#spatial_ref=GEOGCS["WGS 84", DATUM["World Geodetic System 1984", SPHEROID["WGS 84", 6378137.0, 298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich", 0.0, AUTHORITY["EPSG","8901"]], UNIT["degree", 0.017453292519943295], AXIS["Geodetic longitude", EAST], AXIS["Geodetic latitude", NORTH], AUTHORITY["EPSG","4326"]] Discrete_Classification_map#grid_mapping=crs Discrete_Classification_map#long_name=Discrete_Classification_map Discrete_Classification_map#units= Discrete_Classification_map#_FillValue=0 Discrete_Classification_map#_Unsigned=true NC_GLOBAL#Conventions=CF-1.9 NC_GLOBAL#institution=openEO platform NETCDF_DIM_EXTRA={t} NETCDF_DIM_t_DEF={5,4} NETCDF_DIM_t_VALUES={9131,9496,9862,10227,10592} t#axis=T t#long_name=t t#standard_name=t t#units=days since 1990-01-01 x#long_name=longitude x#standard_name=longitude x#units=degrees_east y#long_name=latitude y#standard_name=latitude y#units=degrees_north Corner Coordinates: Upper Left ( 32.3700000, -18.6690000) ( 32d22'12.00"E, 18d40' 8.40"S) Lower Left ( 32.3700000, -20.9412222) ( 32d22'12.00"E, 20d56'28.40"S) Upper Right ( 34.7588889, -18.6690000) ( 34d45'32.00"E, 18d40' 8.40"S) Lower Right ( 34.7588889, -20.9412222) ( 34d45'32.00"E, 20d56'28.40"S) Center ( 33.5644444, -19.8051111) ( 33d33'52.00"E, 19d48'18.40"S) Band 1 Block=256x256 Type=Int16, ColorInterp=Undefined NoData Value=0 Metadata: grid_mapping=crs long_name=Discrete_Classification_map NETCDF_DIM_t=9131 NETCDF_VARNAME=Discrete_Classification_map units= _FillValue=0 _Unsigned=true Band 2 Block=256x256 Type=Int16, ColorInterp=Undefined NoData Value=0 Metadata: grid_mapping=crs long_name=Discrete_Classification_map NETCDF_DIM_t=9496 NETCDF_VARNAME=Discrete_Classification_map units= _FillValue=0 _Unsigned=true Band 3 Block=256x256 Type=Int16, ColorInterp=Undefined NoData Value=0 Metadata: grid_mapping=crs long_name=Discrete_Classification_map NETCDF_DIM_t=9862 NETCDF_VARNAME=Discrete_Classification_map units= _FillValue=0 _Unsigned=true Band 4 Block=256x256 Type=Int16, ColorInterp=Undefined NoData Value=0 Metadata: grid_mapping=crs long_name=Discrete_Classification_map NETCDF_DIM_t=10227 NETCDF_VARNAME=Discrete_Classification_map units= _FillValue=0 _Unsigned=true Band 5 Block=256x256 Type=Int16, ColorInterp=Undefined NoData Value=0 Metadata: grid_mapping=crs long_name=Discrete_Classification_map NETCDF_DIM_t=10592 NETCDF_VARNAME=Discrete_Classification_map units= _FillValue=0 _Unsigned=true
data_array = dataset.ReadAsArray().astype(np.int16)
print(data_array.dtype)
int16
print(data_array.min(), data_array.max())
0 200
plt.figure(figsize=(8,10))
plt.imshow(data_array[0], cmap = 'jet', interpolation='nearest')
plt.title('Land Cover CLasses')
plt.colorbar(shrink=0.4)
plt.grid(False)
plt.show()
#histogram
plt.hist(data_array[0], bins=10)
plt.show()