diff --git a/sen2chain/indices.py b/sen2chain/indices.py
index e415ad6ce00aa914fe50cca3224f57ecc8a8bb1d..624d80ebcc30902b9ddb1c7464273bdd2f98a46b 100644
--- a/sen2chain/indices.py
+++ b/sen2chain/indices.py
@@ -14,7 +14,7 @@ import os
from typing import Union, List
from .indices_functions import (create_raw_ndvi, create_raw_ndwi,
- create_masked_indice, ndvi_tiff_2_jp2)
+ create_masked_indice, index_tiff_2_jp2)
from .colormap import matplotlib_colormap_to_rgb, create_colormap
logger = logging.getLogger(__name__)
@@ -102,7 +102,7 @@ class Ndvi(Indice):
create_raw_ndvi(nir_path=self.l2a_product.b08_10m,
vir_path=self.l2a_product.b04_10m,
out_path=(out_path / self.indice_raw))
- ndvi_tiff_2_jp2(img_path=(out_path / self.indice_raw),
+ index_tiff_2_jp2(img_path=(out_path / self.indice_raw),
out_path=(out_path / self.indice_filename))
if nodata_clouds:
if self.l2a_product.user_cloud_mask is None:
@@ -120,7 +120,7 @@ class Ndvi(Indice):
create_masked_indice(indice_path=ndvi_name,
cloud_mask_path=self.l2a_product.user_cloud_mask,
out_path=(out_path / masked_indice_raw))
- ndvi_tiff_2_jp2(img_path=(out_path / masked_indice_raw),
+ index_tiff_2_jp2(img_path=(out_path / masked_indice_raw),
out_path=(out_path / masked_indice_filename))
os.remove(str(out_path / masked_indice_raw))
@@ -148,7 +148,8 @@ class NdwiGao(Indice):
"""
name = "NDWI_GAO"
filename_template = "{product_identifier}_NDWIg{ext}"
- ext = ".tif"
+ ext = ".jp2"
+ ext_raw = ".tif"
colormap = cm.RdYlBu
# colormap = get_tsv_colormap("../../share/data/RdYlGn.lut")
@@ -162,37 +163,49 @@ class NdwiGao(Indice):
self.out_path = None
# filenames
- self.indice_stem = self.filename_template.format(product_identifier=self.l2a_product.identifier,
- ext="")
+ self.indice_stem = self.filename_template.format(product_identifier=self.l2a_product.identifier, ext="")
self.indice_filename = self.indice_stem + self.ext
+ self.indice_raw = self.indice_stem + self.ext_raw
def process_indice(
self,
out_path: pathlib.PosixPath,
- reprocess: bool = False, nodata_clouds: bool = False, quicklook: bool = False
+ reprocess: bool = False,
+ nodata_clouds: bool = False,
+ quicklook: bool = False
) -> None:
""" process """
self.out_path = out_path
- if (self.out_path / self.indice_filename).exists() and not reprocess:
+ if (out_path / self.indice_filename).exists() and not reprocess:
logger.info("{} already exists".format(self.indice_filename))
else:
create_raw_ndwi(nir_path=self.l2a_product.b08_10m,
swir_path=self.l2a_product.b11_20m,
+ out_path=(out_path / self.indice_raw))
+ index_tiff_2_jp2(img_path=(out_path / self.indice_raw),
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
+ masked_indice_raw = self.indice_stem + "_MASKED" + self.ext_raw
- if (self.out_path / masked_indice_filename).exists() and not reprocess:
+ if (out_path / masked_indice_filename).exists() and not reprocess:
logger.info("{} already exists".format(masked_indice_filename))
else:
- create_masked_indice(indice_path=(self.out_path / self.indice_filename),
+ if (out_path / self.indice_raw).exists():
+ ndwi_name = (out_path / self.indice_raw)
+ else:
+ ndwi_name = (out_path / self.indice_filename)
+ create_masked_indice(indice_path=ndwi_name,
cloud_mask_path=self.l2a_product.user_cloud_mask,
- out_path=(self.out_path / masked_indice_filename))
-
+ out_path=(out_path / masked_indice_raw))
+ index_tiff_2_jp2(img_path=(out_path / masked_indice_raw),
+ out_path=(out_path / masked_indice_filename))
+ os.remove(str(out_path / masked_indice_raw))
+
if quicklook:
cmap = matplotlib_colormap_to_rgb(self.colormap, revers=False)
@@ -200,9 +213,10 @@ class NdwiGao(Indice):
if (self.out_path / quicklook_filename).exists() and not reprocess:
logger.info("{} already exists".format(quicklook_filename))
else:
+ logger.info("creating quicklook")
create_colormap(raster=(self.out_path / self.indice_filename),
cloud_mask=self.l2a_product.user_cloud_mask,
- lut_dict=cmap,
+ lut_dict=cmap, clouds_color="white",
out_path=(self.out_path / quicklook_filename))
diff --git a/sen2chain/indices_functions.py b/sen2chain/indices_functions.py
index 8a4aa5c58e81582011bb28be2ff05a64d6caa334..d59fe8539253e0e23a7c203d4fbbe26b6cc18a26 100644
--- a/sen2chain/indices_functions.py
+++ b/sen2chain/indices_functions.py
@@ -57,13 +57,13 @@ def create_raw_ndwi(nir_path: Union[str, pathlib.PosixPath],
swir_path: Union[str, pathlib.PosixPath],
out_path: Union[str, pathlib.PosixPath]="./raw_ndwi.tif") -> pathlib.PosixPath:
"""
- Creates a NDWI raster from NIR and VIR rasters.
+ Creates a NDWI raster from NIR and SWIR rasters.
:param nir_path: path to the NIR raster.
:param swir_path: path to the SWIR raster.
:param out_path: path to the output raster.
"""
- logger.info("creating raw NDWI (int16)")
+ logger.info("creating raw NDWI (tiff - int16)")
with rasterio.open(str(nir_path)) as nir_src, \
rasterio.open(str(swir_path)) as swir_src:
@@ -78,9 +78,9 @@ def create_raw_ndwi(nir_path: Union[str, pathlib.PosixPath],
swir_reproj = np.empty(nir.shape, dtype=np.float32)
reproject(source=swir,
destination=swir_reproj,
- src_transform=swir_src.affine,
+ src_transform=swir_src.transform,
src_crs=swir_src.crs,
- dst_transform=nir_src.affine,
+ dst_transform=nir_src.transform,
dst_crs=nir_src.crs,
resampling=Resampling.nearest)
@@ -94,11 +94,11 @@ def create_raw_ndwi(nir_path: Union[str, pathlib.PosixPath],
tiled=False,
dtype=np.int16,
nodata=32767,
- transform=nir_src.affine)
-
- with rasterio.open(str(out_path), "w", **nir_profile) as dst:
- dst.write(ndwi_masked, 1)
-
+ transform=nir_src.transform)
+ nir_profile.pop('tiled', None)
+ with rasterio.Env(GDAL_CACHEMAX=512) as env:
+ with rasterio.open(str(out_path), "w", **nir_profile) as dst:
+ dst.write(ndwi_masked, 1)
return Path(str(out_path)).absolute
@@ -150,16 +150,13 @@ def create_masked_indice(
return str(Path(str(out_path)).absolute)
-def ndvi_tiff_2_jp2(img_path: Union[str, pathlib.PosixPath],
- out_path: Union[str, pathlib.PosixPath]="./ndvi_2_jp2.jp2") -> pathlib.PosixPath:
+def index_tiff_2_jp2(img_path: Union[str, pathlib.PosixPath],
+ out_path: Union[str, pathlib.PosixPath]="./indice_2_jp2.jp2") -> pathlib.PosixPath:
"""
- Convert a NDVI file from TIF to JP2.
-
- :param nir_path: path to the NIR raster.
- :param vir_path: path to the VIR raster.
+ Convert a indice file from TIF to JP2.
:param out_path: path to the output raster.
"""
- logger.info("converting raw NDWI to JPEG2000")
+ logger.info("converting raw indice to JPEG2000")
src_ds = gdal.Open(str(img_path))
driver = gdal.GetDriverByName("JP2OpenJPEG")
dst_ds = driver.CreateCopy(str(out_path), src_ds, options=['NBITS=15', 'CODEC=JP2', 'QUALITY=20'])