Skip to content
Snippets Groups Projects
Commit 798063f1 authored by Neha Hunka's avatar Neha Hunka
Browse files

Code for Nature Scientific Data

parent 07966f6c
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Hunka_2024_NSD/NOTES_data_download_and_preprocessing.ipynb
+ 2
− 2
View file @ 798063f1
%% Cell type:markdown id:cf89fb5f-1653-4ec5-aa29-ed7c289890f4 tags:
###################################
The following are notes that provide the steps followed for the download and pre-processing of input datasets (Table 1). The end product are the datasets resampled and aligned to a 30 m x 30 m grid in 10 x 10 degree tiles. These steps are reproducible for batch processing on the AWS DPS cloud-computing system that supports the NASA MAAP. For ease of use, a single 10 x 10 degree tile (00N_000E) is exemplified so the steps are comprehensible and easily implementable on local machines. Cells are executable in either Python or R depending on the opertation to be performed, and it is indicated which to use where.
###################################
%% Cell type:markdown id:2aa16388-c22c-4a4e-8769-06839a3d73d7 tags:
### (1) GLOBAL FOREST HEIGHT LAYERS
#### Downloaded on https://storage.googleapis.com/earthenginepartners-hansen/GFC-2022-v1.10/download.html
#### No pre-processing required
%% Cell type:markdown id:69df931b-998d-4f5a-9dc1-78a6a5da8ce1 tags:
### (2) GLOBAL TREE COVER LOSS YEAR
#### Downloaded on https://storage.googleapis.com/earthenginepartners-hansen/GFC-2022-v1.10/download.html
#### No pre-processing required
%% Cell type:markdown id:66fdfca0-95cf-4db3-8af9-460b1ce4e7c3 tags:
### (3) SPATIAL DATABASE OF PLANTED TREES
#### Download restricted till version 2.0 is officially released (code to be updated when released)
#### No pre-processing required
%% Cell type:markdown id:0faa7be4-854e-41b9-a2f3-93b87880060b tags:
### (4) FOREST LANDSCAPE INTEGRITY INDEX
#### Available at https://www.forestintegrity.com/
%% Cell type:code id:8fdbccf2-9fb1-4e92-b324-973234c26412 tags:
``` python
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
#### DOWNLOADED USING WGET COMMANDS ON TERMINAL, PROCESSED USING GDAL COMMANDS ########
#### PREPROCESSING PROVIDED HERE FOR A SINGLE 10X10 DEGREE TILE ###
TILE <- c("00N_000E")
EXTENT <- c(0,-10,10,0)
gdalUtils::gdalwarp('./Forest_Integrity_Index/FII_world.tif',paste0('./Forest_Integrity_Index/',TILE,'_FII.tif'),tr=c(0.00025,0.00025),a_nodata=0.0,te=EXTENT,r="average",ot="Int16",co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"))
```
%% Cell type:markdown id:8ec571c1-bd24-472b-9ef4-b938f44c899e tags:
### (5) INTACT FORESTS LANDSCAPES
#### Available at https://intactforests.org/
%% Cell type:code id:2784594b-3122-45d0-b1ff-0d0da240dc90 tags:
``` python
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
#### DOWNLOADED USING WGET COMMANDS ON TERMINAL, PROCESSED USING GDAL COMMANDS ########
#### PREPROCESSING PROVIDED HERE FOR A SINGLE 10X10 DEGREE TILE ###
TILE <- c("00N_000E")
EXTENT <- c(0,-10,10,0)
OUT <- paste0('./IFL_rasters_2020/',TILE,'_IFL_2020.tif')
gdalUtils::gdal_rasterize('./IFL_2020.gpkg',OUT,l="ifl_2020",burn=1,tr=c(0.00025,0.00025),a_nodata=0.0,te=EXTENT,ot="Int16",co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"),verbose=TRUE,output_Raster=TRUE)
```
%% Cell type:markdown id:b72fbd80-f2db-4cf0-8d4b-48a2fd2df853 tags:
### (6) PRIMARY HUMID TROPICAL FORESTS
#### Available at https://glad.umd.edu/dataset/primary-forest-humid-tropics
#### Note, data is available for download for each continent covering the tropics. 10 x 10 degree tiles are hence produced separately for Asia, Africa and South America.
%% Cell type:code id:2e7c53d8-ac41-4d96-9e11-f1cc08ab7727 tags:
``` python
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
#### DOWNLOADED USING WGET COMMANDS ON TERMINAL, PROCESSED USING GDAL COMMANDS ########
#### PREPROCESSING PROVIDED HERE FOR A SINGLE 10X10 DEGREE TILE ###
TILE <- c("00N_000E")
EXTENT <- c(0,-10,10,0)
gdalUtils::gdalwarp('./Asia_2001_primary.tif',paste0('./',TILE,'_Asia.tif',tr=c(0.00025,0.00025),a_nodata=0.0,te=EXTENT,ot="Int16",co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"))
```
%% Cell type:markdown id:dddcb298-466b-4f84-8226-88df68c2e48d tags:
### (7) JRC TMF TRANSITION MAP, DEFORESTATION MAP AND DEGRADATION MAP ###
#### Available at https://forobs.jrc.ec.europa.eu/TMF/data#downloads
%% Cell type:code id:d2608083-3051-4392-b792-bb92bd899bf9 tags:
``` python
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
#### DOWNLOADED USING CURL COMMANDS ON TERMINAL, PROCESSED USING GDAL COMMANDS ########
#### Curl example provided for tile lat=N0&lon=E00 #####
# curl -o ./JRC_TMST_00N_000E.tif -L "https://ies-ows.jrc.ec.europa.eu/iforce/tmf_v1/download.py?type=tile&dataset=TransitionMap_Subtypes&lat=N0&lon=E00"
#### PREPROCESSING PROVIDED HERE FOR A SINGLE 10X10 DEGREE TILE ###
TILE <- c("00N_000E")
EXTENT <- c(0,-10,10,0)
gdalUtils::gdalwarp('./JRC_TMST_',TILE,'.tif','./',TILE,'_JRC_Transition_Map.tif',tr=c(0.00025,0.00025),a_nodata=0.0,te=EXTENT,r="near",co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"))
```
%% Cell type:markdown id:db5295d9-5608-43d4-8cd3-d4616f5d6ed5 tags:
### (8) GLOBAL FOREST MANAGEMENT
#### Available at https://pure.iiasa.ac.at/id/eprint/17846/
%% Cell type:code id:c1d36177-63a4-4d6b-b5fc-fca613cfe870 tags:
``` python
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
#### DOWNLOADED USING CURL COMMANDS ON TERMINAL, PROCESSED USING GDAL COMMANDS ########
#### Curl example provided #####
# curl -o ./IIASA/GFM.tif -L "https://zenodo.org/records/4541513/files/FML_v3.2.tif?download=1"
#### PREPROCESSING PROVIDED HERE FOR A SINGLE 10X10 DEGREE TILE ###
TILE <- c("00N_000E")
EXTENT <- c(0,-10,10,0)
gdalUtils::gdalwarp('./IIASA/GFM.tif',paste0('./IIASA/',TILE,'_GFM.tif'),tr=c(0.00025,0.00025),a_nodata=0.0,te=EXTENT,r="near",co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"))
```
%% Cell type:markdown id:f00fc854-44b0-4d49-9eff-f185d071970d tags:
### (9) BOREAL FOREST AGE MAPS
#### Available at https://daac.ornl.gov/ABOVE/guides/Boreal_CanopyCover_StandAge.html
%% Cell type:code id:1d61d2f9-fa9e-40c1-badd-72fac9d51231 tags:
``` python
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
#### DOWNLOADED USING WGET COMMANDS ON TERMINAL, PROCESSED USING GDAL COMMANDS ########
files <- Sys.glob(file.path("./Boreal_CanopyCover_StandAge/data", "*.tif"))
DEST <- "./Boreal_CanopyCover_StandAge/data_INT16/"
for (file in files) {
BASE <- basename(file)
OUT <- paste0(DEST,BASE)
gdalUtils::gdal_translate(file,OUT,ot="Int16",of='GTiff',co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"),verbose=TRUE,output_Raster=TRUE)
}
DEST <- "./Boreal_CanopyCover_StandAge/data_INT16/"
tiles <- Sys.glob(file.path(DEST, "*.tif"))
gdalUtils::mosaic_rasters(gdalfile=tiles,dst_dataset=file.path(DEST,"BOREAL_FOREST_AGES_mosaic.tif"),of="GTiff", gdalwarp_params = list(r = "average",ot="Int16"), co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"), overwrite=TRUE, VERBOSE=TRUE)
#### PREPROCESSING PROVIDED HERE FOR A SINGLE 10X10 DEGREE TILE ###
TILE <- c("00N_000E")
EXTENT <- c(0,-10,10,0)
gdalUtils::gdalwarp('./data_INT16/BOREAL_FOREST_AGES_mosaic.tif',paste0('./BOREAL_AGE/',TILE,'_BOREAL_AGE.tif',tr=c(0.00025,0.00025),a_nodata=0.0,te=EXTENT,r="average",ot="Int16",co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"))
```
%% Cell type:markdown id:7a434b4a-e2ea-49a2-a24f-f8fc51c44d7d tags:
### (10) GEZ BOUNDARY LAYERS
#### GEZ 2010 layer was accessed at: https://www.fao.org/forest-resources-assessment/remote-sensing/global-ecological-zones-gez-mapping/en/
#### Desert ecoregions are removed, the vectors are merged with cotinent boundaries and the output "EcoCont_clean.shp" prodocued in the study is accessible at https://drive.google.com/drive/folders/1l3mrSmR4fFMRcUuKLsOTBDbuFdzYIUfk?usp=drive_link
%% Cell type:code id:ad57d62c-6133-4aa0-9608-626498ed7b3b tags:
``` python
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
#### Downloaded vectors, merged, and rasterized to 30 m x 30 m tiles with a unique burn-in value for each GEZ + Continent combination
#### PREPROCESSING PROVIDED HERE FOR A SINGLE 10X10 DEGREE TILE ###
TILE <- c("00N_000E")
EXTENT <- c(0,-10,10,0)
gdalUtils::gdal_rasterize('./EcoCont_clean.shp',paste0('./EcoCont_v2/',TILE,'_EcoCont.tif'),l="EcoCont_clean",a="EcoCont",tr=c(0.00025,0.00025),a_nodata=0.0,te=EXTENT,ot="Int16",co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"),verbose=TRUE,output_Raster=TRUE)
```
%% Cell type:markdown id:ec04cd3f-2279-4e44-b7ca-651cfa727804 tags:
### (11) ESA WORLD COVER 2021 DATASET
#### The ESA WorldCover 2021 dataset was accessed through the NASA MAAP STAC end-point. Steps for this are provided below, and executed in Python. Users outside the NASA MAAP are directed to donwload the original data at https://worldcover2021.esa.int/download, and then follow the pre-processing steps described in the next cell.
%% Cell type:code id:95bd5bd5-a8e4-43dd-8f82-5e20b9b1bbaf tags:
``` python
################ Python 3.10.13 | packaged by conda-forge ################
################ Python 3.10.13 | packaged by conda-forge ################
################ Python 3.10.13 | packaged by conda-forge ################
########## ACCESS AND DOWNLOAD THROUGH MAAP STAC ##########
import glob
import stackstac
from zipfile import ZipFile
from stackstac import stack, mosaic
import pystac_client
import os
from rio_tiler.io import STACReader
import rioxarray
from rioxarray import merge
import matplotlib.pyplot as plt
os.environ['AWS_NO_SIGN_REQUEST'] = 'YES'
gdal_env=stackstac.DEFAULT_GDAL_ENV.updated({'AWS_NO_SIGN_REQUEST': 'YES'})
URL = "https://services.terrascope.be/stac/"
catalog = pystac_client.Client.open(URL)
xmin = -180
ymin = -70
step = 3
nx = 122
ny = 50
for x in range(0,nx-1):
for y in range(0,ny-1):
IDEN = str((xmin+(step*x))+1.4)+"_"+str((ymin+(step*y))+1.4)
if not os.path.exists("./ESA_WorldCover2021_v2/WC_"+IDEN+".tif"):
bbox = [(xmin+(step*x))+1.4,(ymin+(step*y))+1.4,(xmin+(step*(x+1)))-1.4,(ymin+(step*(y+1)))-1.4]
stac_collection = catalog.search(
collections=["urn:eop:VITO:ESA_WorldCover_10m_2021_AWS_V2"],
bbox=bbox
)
try:
arr = stack(stac_collection.get_all_items(),epsg=4326,resolution=0.00008333333333333333055,gdal_env=gdal_env)
arr_est = arr.sel(band='ESA_WORLDCOVER_10M_MAP')
arr_masked_lt = arr_est.where(arr_est == 10)
arr_masked_lt[0].rio.to_raster(raster_path="./ESA_WorldCover2021_v2/WC_"+IDEN+".tif", driver='GTiff', compress='lzw',dtype="int8")
except:
print(IDEN)
```
%% Cell type:code id:b90c6063-b077-499e-8bfd-b555b21f4421 tags:
``` python
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
################ R version 4.2.3 ###########################
DEST <- "./ESA_WorldCover2021_v2/"
tiles <- Sys.glob(file.path(DEST, "WC_*.tif"))
tiles <- tiles[!str_detect(tiles,pattern="30m.tif")]
tiles <- tiles[!str_detect(tiles,pattern="S_")]
tiles <- tiles[!str_detect(tiles,pattern="N_")]
tiles <- tiles[!str_detect(tiles,pattern="E_")]
tiles <- tiles[!str_detect(tiles,pattern="W_")]
tiles <- tiles[!str_detect(tiles,pattern="W.tif")]
tiles <- tiles[!str_detect(tiles,pattern="E.tif")]
for (tile in tiles) {
out_file <- file.path(DEST,paste0(tools::file_path_sans_ext(basename(tile)),"_30m.tif"))
if (!file.exists(out_file)) {
gdalUtils::gdalwarp(tile,out_file,tr=c(0.00025,0.00025),a_nodata=0.0,r="mode",ot="Byte",co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"))
}
}
packages <- c("raster","stringr", "rgdal", "terra", "rgeos", "geosphere","RColorBrewer","gdalUtils","parallel","snow","plyr","maptools", "scales", "sf", "sp","foreign","aws.s3","aws.ec2metadata")
# packages <- c("raster", "rgdal", "terra", "sf", "sp")
package.check <- lapply(packages, FUN = function(x) {
if (!require(x, character.only = TRUE)) {
# install.packages(x, dependencies = TRUE)
library(x, character.only = TRUE, quietly=TRUE)
library(geojsonio)
}
})
DEST <- "./ESA_WorldCover2021_v2/"
xseq <- seq(from = 90, to = 100, by = 10)
yseq <- seq(from = -10, to = 0, by = 10)
for (n in 1:(length(xseq)-1)){
for (m in 1:(length(yseq)-1)){
tiles <- Sys.glob(file.path(DEST, "WC_*_30m.tif"))
tile_xy = c()
for (tile in tiles) {
tile_x <- as.double(strsplit(basename(tile), split = "_")[[1]][2])
tile_y <- as.double(strsplit(basename(tile), split = "_")[[1]][3])
if ((tile_x >= (xseq[n]-4)) & (tile_x < (xseq[n+1]+4)) & (tile_y >= (yseq[m]-10)) & (tile_y < yseq[m+1])) {tile_xy <-append(tile_xy, tile)}
}
if (!is.null(tile_xy)) {
if (yseq[m+1]>=0) {NorS <- "N"}
if (yseq[m+1]<0) {NorS <- "S"}
if (xseq[n]<0) {EorW <- "W"}
if (xseq[n]>=0) {EorW <- "E"}
if (nchar(trunc(abs(xseq[n])))==1){EorWvalue <- paste0("00",as.character(abs(xseq[n])))}
if (nchar(trunc(abs(xseq[n])))==2){EorWvalue <- paste0("0",as.character(abs(xseq[n])))}
if (nchar(trunc(abs(xseq[n])))==3){EorWvalue <- paste0(as.character(abs(xseq[n])))}
if (nchar(trunc(abs(yseq[m+1])))==1){NorSvalue <- paste0("0",as.character(yseq[m+1]))}
if (nchar(trunc(abs(yseq[m+1])))==2){NorSvalue <- paste0(as.character(yseq[m+1]))}
OUT_STRING <- paste0(as.character(NorSvalue),NorS,"_",EorWvalue,EorW)
if (!file.exists(file.path(DEST,paste0("WC_",OUT_STRING,"_C.tif")))){
gdalUtils::mosaic_rasters(gdalfile=tile_xy,dst_dataset=file.path(DEST,paste0("WC_",OUT_STRING,".tif")),of="GTiff", gdalwarp_params = list(r = "max",ot="Byte"), co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"), overwrite=TRUE, VERBOSE=TRUE)
gdalUtils::gdalwarp(file.path(DEST,paste0("WC_",OUT_STRING,".tif")),file.path(DEST,paste0("WC_",OUT_STRING,"_C.tif")),tr=c(0.00025,0.00025),a_nodata=0.0,te=c(xseq[n],yseq[m],xseq[n+1],yseq[m+1]),r="average",ot="Int16",co=c("COMPRESS=DEFLATE","PREDICTOR=2","ZLEVEL=9"))
} else {print("done")}
}
}
}
```
%% Cell type:markdown id:b859af42-1794-4145-9315-a9d084576a69 tags:
### HOW TO HANDLE MISSING TILES
#### All sourced input datasets might not contain the every 10 x 10 degree tile to cover the entire globe. For example, the Primary Humid Tropical Forest layer is provided only for the tropics. As such, for missing tiles, psuedo layers with 0 values are created to cover the entire land surface of the earth. These pseudo layers ensure that the forest classification can run without returning errors for any 10 x 10 degree tile in the world.
#### All sourced input datasets might not contain every 10 x 10 degree tile to cover the entire globe. For example, the Primary Humid Tropical Forest layer is provided only for the tropics. As such, for missing tiles, psuedo layers with 0 values are created to cover the entire land surface of the earth. These pseudo layers ensure that the forest classification can run without returning errors for any 10 x 10 degree tile in the world.
#### Below is an example of such layers created for where the JRC Tropcal maps have no detected degradation
#### Below is an example of such psuedo layers created for where the JRC Tropical maps that are missing from our global-coverage
%% Cell type:code id:a0998db3-d1eb-4e30-bb1f-678fd3f69c08 tags:
``` python
################ Python 3.10.13 | packaged by conda-forge ################
################ Python 3.10.13 | packaged by conda-forge ################
################ Python 3.10.13 | packaged by conda-forge ################
####### CREATE EMPTY JRC LAYERS (only 0 values) FOR MISSING TILES ###########
Blank = rasterio.open("./JRC_DEG/00N_000E_JRC_Degradation_Map.tif").read(1)
Blank = Blank*0
kwargs = rasterio.open("./JRC_DEG/00N_000E_JRC_Degradation_Map.tif").meta
kwargs.update(dtype=rasterio.int16,count=1,compress='lzw')
TILES = ["00N_000E","00N_010E","00N_020E","00N_030E","00N_040E","00N_040W","00N_050E","00N_050W","00N_060W","00N_070E","00N_070W","00N_080W","00N_090E","00N_090W","00N_100E","00N_100W","00N_110E","00N_120E","00N_130E","00N_140E","00N_150E","00N_160E","00N_170E","10N_000E","10N_010E","10N_010W","10N_020E","10N_020W","10N_030E","10N_040E","10N_050E","10N_050W","10N_060W","10N_070E","10N_070W","10N_080E","10N_080W","10N_090E","10N_090W","10N_100E","10N_100W","10N_110E","10N_120E","10N_130E","10N_150E","10N_160E","10N_170E","10S_010E","10S_020E","10S_030E","10S_040E","10S_040W","10S_050E","10S_050W","10S_060W","10S_070W","10S_080W","10S_110E","10S_120E","10S_130E","10S_140E","10S_150E","10S_160E","10S_170E","10S_180W","20N_000E","20N_010E","20N_010W","20N_020E","20N_020W","20N_030E","20N_030W","20N_040E","20N_050E","20N_060W","20N_070E","20N_070W","20N_080E","20N_080W","20N_090E","20N_090W","20N_100E","20N_100W","20N_110E","20N_110W","20N_120E","20N_120W","20N_160W","20S_010E","20S_020E","20S_030E","20S_040E","20S_050E","20S_050W","20S_060W","20S_070W","20S_080W","20S_110E","20S_120E","20S_130E","20S_140E","20S_150E","20S_160E","20S_180W","30N_000E","30N_010E","30N_010W","30N_020E","30N_020W","30N_030E","30N_040E","30N_050E","30N_060E","30N_070E","30N_080E","30N_080W","30N_090E","30N_090W","30N_100E","30N_100W","30N_110E","30N_110W","30N_120E","30N_120W","30N_130E","30N_160W","30N_170W","30S_010E","30S_020E","30S_030E","30S_060W","30S_070W","30S_080W","30S_110E","30S_120E","30S_130E","30S_140E","30S_150E","30S_170E","40N_000E","40N_010E","40N_010W","40N_020E","40N_020W","40N_030E","40N_040E","40N_050E","40N_060E","40N_070E","40N_070W","40N_080E","40N_080W","40N_090E","40N_090W","40N_100E","40N_100W","40N_110E","40N_110W","40N_120E","40N_120W","40N_130E","40N_130W","40N_140E","40S_070W","40S_080W","40S_140E","40S_160E","40S_170E","50N_000E","50N_010E","50N_010W","50N_020E","50N_030E","50N_040E","50N_050E","50N_060E","50N_060W","50N_070E","50N_070W","50N_080E","50N_080W","50N_090E","50N_090W","50N_100E","50N_100W","50N_110E","50N_110W","50N_120E","50N_120W","50N_130E","50N_130W","50N_140E","50N_150E","50S_060W","50S_070W","50S_080W","60N_000E","60N_010E","60N_010W","60N_020E","60N_020W","60N_030E","60N_040E","60N_050E","60N_060E","60N_060W","60N_070E","60N_070W","60N_080E","60N_080W","60N_090E","60N_090W","60N_100E","60N_100W","60N_110E","60N_110W","60N_120E","60N_120W","60N_130E","60N_130W","60N_140E","60N_140W","60N_150E","60N_150W","60N_160E","60N_160W","60N_170E","60N_170W","60N_180W","70N_000E","70N_010E","70N_020E","70N_030E","70N_040E","70N_050E","70N_060E","70N_070E","70N_070W","70N_080E","70N_080W","70N_090E","70N_090W","70N_100E","70N_100W","70N_110E","70N_110W","70N_120E","70N_120W","70N_130E","70N_130W","70N_140E","70N_140W","70N_150E","70N_150W","70N_160E","70N_160W","70N_170E","70N_170W","70N_180W","80N_010E","80N_020E","80N_030E","80N_070E","80N_080E","80N_090E","80N_100E","80N_110E","80N_120E","80N_130E","80N_130W","80N_140E","80N_140W","80N_150E","80N_150W","80N_160E","80N_160W","80N_170E","80N_170W"]
for TILE in TILES:
if not os.path.exists("./JRC_DEG/"+TILE+"_JRC_Degradation_Map.tif"):
with rasterio.open("./JRC_DEG/"+TILE+"_JRC_Degradation_Map.tif", 'w', **kwargs) as dst: dst.write_band(1, Blank.astype(rasterio.int16))
```
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment