Skip to content

Commit

Permalink
Merge branch 'Johannes-R-Schmid-development' into review-pr-34
Browse files Browse the repository at this point in the history
  • Loading branch information
@devisperessutti
devisperessutti committed Jan 6, 2019
2 parents ad8af3b + 2879e78 commit 3878997
Show file tree
Hide file tree
Showing 3 changed files with 161 additions and 13 deletions.
Binary file added example_data/TestEOPatch/data/REFERENCE_SCENES.npy
View file
Binary file not shown.
27 changes: 14 additions & 13 deletions features/eolearn/features/radiometric_normalization.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -123,27 +123,27 @@ def _geoville_index_by_percentile(self, data, percentile):
def _get_reference_band(self, data):
""" Extract reference band from input 4D data according to compositing method
:param data: 4D array from which to extract reference band (e.g. blue, maxNDVI, ..)
:type data: numpy array
:return: 3D array containing reference band according to compositing method
:param data: 4D array from which to extract reference band (e.g. blue, maxNDVI, ..)
:type data: numpy array
:return: 3D array containing reference band according to compositing method
"""
raise NotImplementedError

def _get_indices(self, data):
""" Compute indices along temporal dimension corresponding to the sought percentile
:param data: Input 3D array holding the reference band
:type data: numpy array
:return: 2D array holding the temporal index corresponding to percentile
:param data: Input 3D array holding the reference band
:type data: numpy array
:return: 2D array holding the temporal index corresponding to percentile
"""
indices = self._index_by_percentile(data, self.percentile)
return indices

def execute(self, eopatch):
""" Compute composite array merging temporal frames according to the compositing method
:param eopatch: eopatch holding time-series
:return: eopatch with composite image of time-series
:param eopatch: eopatch holding time-series
:return: eopatch with composite image of time-series
"""
feature_type, feature_name = next(self.feature(eopatch))
data = eopatch[feature_type][feature_name].copy()
Expand Down Expand Up @@ -310,9 +310,9 @@ def _get_reference_band(self, data):
The max-ratio is defined as max(NIR,SWIR1)/BLUE
:param data: 4D array from which to compute the max-ratio reference band
:type data: numpy array
:return: 3D array containing the max-ratio reference band
:param data: 4D array from which to compute the max-ratio reference band
:type data: numpy array
:return: 3D array containing the max-ratio reference band
"""
blue = data[..., self.blue_idx].astype("float32")
nir = data[..., self.nir_idx].astype("float32")
Expand Down Expand Up @@ -343,8 +343,9 @@ def __init__(self, feature, reference):
def execute(self, eopatch):
""" Perform histogram matching of the time-series with respect to a reference scene
:param eopatch: eopatch holding the time-series and reference data
:return: The same eopatch instance with the normalised time-series
:param eopatch: eopatch holding the time-series and reference data
:type eopatch: EOPatch
:return: The same eopatch instance with the normalised time-series
"""
feature_type, feature_name, new_feature_name = next(self.feature(eopatch))
reference_type, reference_name = next(self.reference(eopatch))
Expand Down
147 changes: 147 additions & 0 deletions features/eolearn/tests/test_radiometric_normalization.py
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,147 @@
import unittest
import os.path
import numpy as np
from datetime import datetime

from eolearn.core import EOPatch, FeatureType
from eolearn.mask import MaskFeature
from eolearn.features import ReferenceScenes, BlueCompositing, HOTCompositing, MaxNDVICompositing, \
MaxNDWICompositing, MaxRatioCompositing, HistogramMatching


class TestRadiometricNormalization(unittest.TestCase):

class RadiometricNormalizationTestCase:
"""
Container for each interpolation test case
"""
TEST_PATCH_FILENAME = os.path.join(os.path.dirname(os.path.realpath(__file__)), '../../../example_data',
'TestEOPatch')

def __init__(self, name, task, img_min=None, img_max=None, img_mean=None, img_median=None):
self.name = name
self.task = task
self.img_min = img_min
self.img_max = img_max
self.img_mean = img_mean
self.img_median = img_median
self.result = None

def execute(self):
patch = EOPatch.load(self.TEST_PATCH_FILENAME)
np.random.seed(0)
patch.mask['SCL'] = np.random.randint(0, 11, patch.data['BANDS-S2-L1C'].shape, np.uint8)
blue = BlueCompositing((FeatureType.DATA, 'REFERENCE_SCENES'),
(FeatureType.DATA_TIMELESS, 'REFERENCE_COMPOSITE'), blue_idx=0,
interpolation='geoville')
blue.execute(patch)
self.result = self.task.execute(patch)

@classmethod
def setUpClass(cls):

cls.test_cases = [
cls.RadiometricNormalizationTestCase('mask feature', MaskFeature(
(FeatureType.DATA, 'BANDS-S2-L1C', 'TEST'), (FeatureType.MASK, 'SCL'),
mask_values=[0, 1, 2, 3, 8, 9, 10, 11]), img_min=0.0002, img_max=1.4244, img_mean=0.21167801,
img_median=0.142),

cls.RadiometricNormalizationTestCase('reference scene', ReferenceScenes(
(FeatureType.DATA, 'BANDS-S2-L1C', 'TEST'), (FeatureType.SCALAR, 'CLOUD_COVERAGE'),
max_scene_number=5), img_min=0.0005, img_max=0.5318, img_mean=0.16823094,
img_median=0.1404),

cls.RadiometricNormalizationTestCase('blue compositing', BlueCompositing(
(FeatureType.DATA, 'REFERENCE_SCENES'), (FeatureType.DATA_TIMELESS, 'TEST'),
blue_idx=0, interpolation='geoville'), img_min=0.0005, img_max=0.5075, img_mean=0.11658352,
img_median=0.0833),

cls.RadiometricNormalizationTestCase('hot compositing', HOTCompositing(
(FeatureType.DATA, 'REFERENCE_SCENES'), (FeatureType.DATA_TIMELESS, 'TEST'),
blue_idx=0, red_idx=2, interpolation='geoville'), img_min=0.0005, img_max=0.5075, img_mean=0.117758796,
img_median=0.0846),

cls.RadiometricNormalizationTestCase('max ndvi compositing', MaxNDVICompositing(
(FeatureType.DATA, 'REFERENCE_SCENES'), (FeatureType.DATA_TIMELESS, 'TEST'),
red_idx=2, nir_idx=7, interpolation='geoville'), img_min=0.0005, img_max=0.5075,
img_mean=0.13430128, img_median=0.0941),

cls.RadiometricNormalizationTestCase('max ndwi compositing', MaxNDWICompositing(
(FeatureType.DATA, 'REFERENCE_SCENES'), (FeatureType.DATA_TIMELESS, 'TEST'),
nir_idx=6, swir1_idx=8, interpolation='geoville'), img_min=0.0005, img_max=0.5318,
img_mean=0.2580135, img_median=0.2888),

cls.RadiometricNormalizationTestCase('max ratio compositing', MaxRatioCompositing(
(FeatureType.DATA, 'REFERENCE_SCENES'), (FeatureType.DATA_TIMELESS, 'TEST'),
blue_idx=0, nir_idx=6, swir1_idx=8, interpolation='geoville'), img_min=0.0006, img_max=0.5075,
img_mean=0.13513365, img_median=0.0958),

cls.RadiometricNormalizationTestCase('histogram matching', HistogramMatching(
(FeatureType.DATA, 'BANDS-S2-L1C', 'TEST'),
(FeatureType.DATA_TIMELESS, 'REFERENCE_COMPOSITE')),
img_min=-0.049050678, img_max=0.68174845, img_mean=0.1165936, img_median=0.08370649)
]

for test_case in cls.test_cases:
test_case.execute()

def test_output_type(self):
for test_case in self.test_cases:
with self.subTest(msg='Test case {}'.format(test_case.name)):
self.assertTrue(isinstance(test_case.result.timestamp, list), "Expected a list of timestamps")
self.assertTrue(isinstance(test_case.result.timestamp[0], datetime),
"Expected timestamps of type datetime.datetime")

def test_stats(self):
for tci, test_case in enumerate(self.test_cases):
delta = 1e-3

if tci not in [2, 3, 4, 5, 6]:
if test_case.img_min is not None:
min_val = np.nanmin(test_case.result.data['TEST'])
with self.subTest(msg='Test case {}'.format(test_case.name)):
self.assertAlmostEqual(test_case.img_min, min_val, delta=delta,
msg="Expected min {}, got {}".format(test_case.img_min, min_val))
if test_case.img_max is not None:
max_val = np.nanmax(test_case.result.data['TEST'])
with self.subTest(msg='Test case {}'.format(test_case.name)):
self.assertAlmostEqual(test_case.img_max, max_val, delta=delta,
msg="Expected max {}, got {}".format(test_case.img_max, max_val))
if test_case.img_mean is not None:
mean_val = np.nanmean(test_case.result.data['TEST'])
with self.subTest(msg='Test case {}'.format(test_case.name)):
self.assertAlmostEqual(test_case.img_mean, mean_val, delta=delta,
msg="Expected mean {}, got {}".format(test_case.img_mean, mean_val))
if test_case.img_median is not None:
median_val = np.nanmedian(test_case.result.data['TEST'])
with self.subTest(msg='Test case {}'.format(test_case.name)):
self.assertAlmostEqual(test_case.img_median, median_val, delta=delta,
msg="Expected median {}, got {}".format(test_case.img_median,
median_val))

else:
if test_case.img_min is not None:
min_val = np.nanmin(test_case.result.data_timeless['TEST'])
with self.subTest(msg='Test case {}'.format(test_case.name)):
self.assertAlmostEqual(test_case.img_min, min_val, delta=delta,
msg="Expected min {}, got {}".format(test_case.img_min, min_val))
if test_case.img_max is not None:
max_val = np.nanmax(test_case.result.data_timeless['TEST'])
with self.subTest(msg='Test case {}'.format(test_case.name)):
self.assertAlmostEqual(test_case.img_max, max_val, delta=delta,
msg="Expected max {}, got {}".format(test_case.img_max, max_val))
if test_case.img_mean is not None:
mean_val = np.nanmean(test_case.result.data_timeless['TEST'])
with self.subTest(msg='Test case {}'.format(test_case.name)):
self.assertAlmostEqual(test_case.img_mean, mean_val, delta=delta,
msg="Expected mean {}, got {}".format(test_case.img_mean, mean_val))
if test_case.img_median is not None:
median_val = np.nanmedian(test_case.result.data_timeless['TEST'])
with self.subTest(msg='Test case {}'.format(test_case.name)):
self.assertAlmostEqual(test_case.img_median, median_val, delta=delta,
msg="Expected median {}, got {}".format(test_case.img_median,
median_val))


if __name__ == '__main__':
unittest.main()

0 comments on commit 3878997

Please sign in to comment.