From 29e39ececd14bbf4930140652cdac63ea1167fa2 Mon Sep 17 00:00:00 2001
From: Impact <pascal.mouquet@ird.fr>
Date: Fri, 30 Nov 2018 16:28:05 +0400
Subject: [PATCH] TIFF / JPEG2000 & NDVI adjustments
---
sen2chain/cloud_mask.py | 28 +++++++++++-----------------
sen2chain/indices.py | 2 --
sen2chain/indices_functions.py | 23 ++++++++++++-----------
sen2chain/products.py | 2 +-
4 files changed, 24 insertions(+), 31 deletions(-)
diff --git a/sen2chain/cloud_mask.py b/sen2chain/cloud_mask.py
index cb5d68f..437cb11 100755
--- a/sen2chain/cloud_mask.py
+++ b/sen2chain/cloud_mask.py
@@ -3,6 +3,7 @@
"""
Module for computing cloud-mask from a sen2cor clouds classification raster.
"""
+
import logging
import pathlib
from pathlib import Path
@@ -16,7 +17,8 @@ from rasterio.features import shapes, rasterize
import numpy as np
from typing import Mapping, Sequence, AbstractSet, List, Set, Dict, Tuple, Optional, Union
import gdal
-
+import time
+import os
# type annotations
from .config import Config
@@ -254,14 +256,10 @@ def create_cloud_mask_v2(cloud_mask: Union[str, pathlib.PosixPath], out_path="./
LOGGER.info("Creating cloud-mask_v2")
out_temp_path = Path(Config().get("temp_path"))
- # LOGGER.info(out_temp_path)
- # LOGGER.info(str(out_temp_path))
- # ATTENTION, si jamais des traitements en parallèle: les fichiers
- # temporaires auront le même nom => conflits
- # Ã utiliser uniquement pour des tests
- out_bin = str(out_temp_path / "temp_bin.tif")
- out_erode = str(out_temp_path / "temp_erode.tif")
- out_dilate = str(out_temp_path / "temp_dilate.tif")
+ out_mil_epoch = str(int(time.time()*1000))
+ out_bin = str(out_temp_path / Path("tmp_bin_" + out_mil_epoch + ".tif"))
+ out_erode = str(out_temp_path / Path("tmp_erode_" + out_mil_epoch + ".tif"))
+ out_dilate = str(out_temp_path / Path("tmp_dilate" + out_mil_epoch + ".tif"))
# Seuillage du masque nuage Sen2Cor
CLD_seuil = 25
@@ -271,7 +269,6 @@ def create_cloud_mask_v2(cloud_mask: Union[str, pathlib.PosixPath], out_path="./
BandMath.SetParameterOutputImagePixelType("out", 0)
BandMath.SetParameterString("exp", "im1b1 > " + str(CLD_seuil))
BandMath.ExecuteAndWriteOutput()
- LOGGER.info("bin fini")
# Erosion
BinaryMorphologicalOperation = otbApplication.Registry.CreateApplication("BinaryMorphologicalOperation")
@@ -283,9 +280,7 @@ def create_cloud_mask_v2(cloud_mask: Union[str, pathlib.PosixPath], out_path="./
BinaryMorphologicalOperation.SetParameterInt("structype.ball.yradius", 1)
BinaryMorphologicalOperation.SetParameterFloat("filter.erode.foreval", 1)
BinaryMorphologicalOperation.SetParameterString("filter","erode")
- LOGGER.info("erosion debut")
BinaryMorphologicalOperation.ExecuteAndWriteOutput()
- LOGGER.info("erosion fini")
# Dilatation
BinaryMorphologicalOperation = otbApplication.Registry.CreateApplication("BinaryMorphologicalOperation")
@@ -298,16 +293,15 @@ def create_cloud_mask_v2(cloud_mask: Union[str, pathlib.PosixPath], out_path="./
BinaryMorphologicalOperation.SetParameterFloat("filter.dilate.foreval", 1)
BinaryMorphologicalOperation.SetParameterString("filter","dilate")
BinaryMorphologicalOperation.ExecuteAndWriteOutput()
- LOGGER.info("dilatation fini")
- # Enregistrement JP2000
+ # Save to JP2000
src_ds = gdal.Open(out_dilate)
driver = gdal.GetDriverByName("JP2OpenJPEG")
dst_ds = driver.CreateCopy(str(out_path), src_ds, options=["CODEC=JP2", "QUALITY=100", "REVERSIBLE=YES", "YCBCR420=NO"])
dst_ds = None
src_ds = None
- #os.remove(out_bin)
- #os.remove(out_erode)
- #os.remove(out_dilate)
+ os.remove(out_bin)
+ os.remove(out_erode)
+ os.remove(out_dilate)
diff --git a/sen2chain/indices.py b/sen2chain/indices.py
index 07c97e4..d86f672 100644
--- a/sen2chain/indices.py
+++ b/sen2chain/indices.py
@@ -86,11 +86,9 @@ class Ndvi(Indice):
out_path=(out_path / self.indice_filename))
if nodata_clouds:
-
if self.l2a_product.user_cloud_mask is None:
raise ValueError("Cloud mask does not exist")
masked_indice_filename = self.indice_stem + "_MASKED" + self.ext
-
if (self.out_path / masked_indice_filename).exists() and not reprocess:
logger.info("{} already exists".format(masked_indice_filename))
else:
diff --git a/sen2chain/indices_functions.py b/sen2chain/indices_functions.py
index 3d2cae5..c9017be 100644
--- a/sen2chain/indices_functions.py
+++ b/sen2chain/indices_functions.py
@@ -29,27 +29,28 @@ def create_raw_ndvi(nir_path: Union[str, pathlib.PosixPath],
LOGGER.info("creating raw NDVI (int16)")
with rasterio.open(str(nir_path)) as nir_src, \
rasterio.open(str(vir_path)) as vir_src:
-
+ #LOGGER.info("001")
nir_profile = nir_src.profile
# vir_profile = vir_src.profile
-
+ #print(nir_profile)
+ #print(nir_src.transform)
np.seterr(divide='ignore', invalid='ignore') # ignore warnings when dividing by zero
nir = nir_src.read(1).astype(np.float32)
vir = vir_src.read(1).astype(np.float32)
ndvi = ((nir - vir) / (nir + vir)*10000).astype(np.int16)
ndvi_masked = np.where(nir != 0, ndvi, 32767)
-
+ #LOGGER.info("002")
nir_profile.update(driver="Gtiff",
compress="DEFLATE",
tiled=False,
dtype=np.int16,
nodata=32767,
- transform=nir_src.affine)
-
+ transform=nir_src.transform)
+ #LOGGER.info("003")
with rasterio.open(str(out_path), "w", **nir_profile) as dst:
dst.write(ndvi_masked, 1)
-
+ #LOGGER.info("004")
return Path(str(out_path)).absolute
@@ -123,22 +124,22 @@ def create_masked_indice(raw_indice_path: Union[str, pathlib.PosixPath],
profile = indice_src.profile
raw_indice = indice_src.read(1)
cld = cld_src.read(1)
-
+ LOGGER.info("001")
# repoject cloud_mask to ndvi resolution
cld_reproj = np.empty(raw_indice.shape, dtype=np.uint8)
reproject(source=cld,
destination=cld_reproj,
- src_transform=cld_src.affine,
+ src_transform=cld_src.transform,
src_crs=cld_src.crs,
- dst_transform=indice_src.affine,
+ dst_transform=indice_src.transform,
dst_crs=indice_src.crs,
resampling=Resampling.nearest)
-
+ LOGGER.info("002")
indice_borders_mask = np.where(raw_indice > 0, raw_indice, 32767)
indice_cloud_mask = np.where(cld_reproj == 0, indice_borders_mask, 32767)
profile.update(nodata=32767,
- transform=indice_src.affine)
+ transform=indice_src.transform)
with rasterio.open(str(out_path), "w", **profile) as dst:
dst.write(indice_cloud_mask, 1)
diff --git a/sen2chain/products.py b/sen2chain/products.py
index 6516d71..409625c 100755
--- a/sen2chain/products.py
+++ b/sen2chain/products.py
@@ -345,7 +345,7 @@ class L2aProduct(Product):
self.indices_path = self._indices_library_path
# user cloud mask
- self.user_cloud_mask = self.path.parent / (self.identifier + "_CLOUD_MASK.tif")
+ self.user_cloud_mask = self.path.parent / (self.identifier + "_CLOUD_MASK.jp2")
if not self.user_cloud_mask.exists():
self.user_cloud_mask = None
--
GitLab