From 0943add3a874128d738080ed9c7f781f577c4a22 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Hollstein?= Date: Wed, 30 Aug 2017 14:11:33 +0200 Subject: [PATCH 1/3] Split vnir and swir models. --- enpt/io/reader.py | 9 ++++++--- enpt/model/metadata.py | 20 +++++++++++++++++++- 2 files changed, 25 insertions(+), 4 deletions(-) diff --git a/enpt/io/reader.py b/enpt/io/reader.py index 7c59814..385af79 100644 --- a/enpt/io/reader.py +++ b/enpt/io/reader.py @@ -17,7 +17,8 @@ class L1B_Reader(object): self.cfg = user_inputs @staticmethod - def read_inputdata(root_dir, observation_time: datetime, lon_lat_smpl=(15, 15), nsmile_coef=5): + def read_inputdata(root_dir, observation_time: datetime, lon_lat_smpl=(15, 15), nsmile_coef=5, + snr_vnir=None, snr_swir=None): # TODO move to init? """Read L1B, DEM and spatial reference data. @@ -50,8 +51,10 @@ class L1B_Reader(object): # compute radiance and calculate snr L1_obj.DN2TOARadiance() - L1_obj.vnir.detector_meta.calc_snr(data=L1_obj.vnir.data) - L1_obj.swir.detector_meta.calc_snr(data=L1_obj.swir.data) + if snr_swir is not None: + L1_obj.vnir.detector_meta.calc_snr_vnir(data=L1_obj.vnir, snr_data_fn=snr_vnir) + if snr_swir is not None: + L1_obj.swir.detector_meta.calc_snr_swir(data=L1_obj.swir, snr_data_fn=snr_swir) return L1_obj diff --git a/enpt/model/metadata.py b/enpt/model/metadata.py index 011602d..3d4dd33 100644 --- a/enpt/model/metadata.py +++ b/enpt/model/metadata.py @@ -6,6 +6,9 @@ from xml.etree import ElementTree import logging import numpy as np from scipy.interpolate import interp2d +import spectral as sp + +from ..model.images import EnMAP_Detector_SensorGeo L1B_product_props = dict( xml_detector_label=dict( @@ -105,6 +108,7 @@ class EnMAP_Metadata_L1B_Detector_SensorGeo(object): self.lons = self.interpolate_corners(*self.lon_UL_UR_LL_LR, *lon_lat_smpl) self.unit = 'none' # '" ".join(xml.findall("%s/radiance_unit" % lbl)[0].text.split()) self.unitcode = 'DN' + self.snr = None def calc_smile(self): """Compute smile for each EnMAP column. @@ -120,7 +124,21 @@ class EnMAP_Metadata_L1B_Detector_SensorGeo(object): self.smile_coef # shape = (nwvl, nsmile_coef) ) # shape = (ncols, nwvl) - def calc_snr(self, data: np.ndarray): + def calc_snr_vnir(self, data: EnMAP_Detector_SensorGeo, snr_data_fn: str): + """Compute EnMAP SNR from radiance data. + + :param data: Numpy array with radiance for scene + """ + self.logger.info("Compute snr for: %s" % self.detector_name) + ds = sp.open_image(snr_data_fn) + p_hg = ds[0:3, :, :] + p_lg = ds[0:3, :, :] + l_th = ds[6, :, :] + + + return 500 * np.ones(data.shape, dtype=np.float) + + def calc_snr_swir(self, data): """Compute EnMAP SNR from radiance data. :param data: Numpy array with radiance for scene -- GitLab From 1647339aa8e7d1994b6c048f96b7a929126ac469 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Hollstein?= Date: Wed, 30 Aug 2017 16:30:43 +0200 Subject: [PATCH 2/3] Add snr model to input reader. --- enpt/io/reader.py | 11 +++--- enpt/model/metadata.py | 35 ++++++++++++------- .../data/EnMAP_Sensor/EnMAP_Level_1B_SNR.zip | 3 ++ tests/test_l1b_reader.py | 21 ++++++++--- 4 files changed, 47 insertions(+), 23 deletions(-) create mode 100644 tests/data/EnMAP_Sensor/EnMAP_Level_1B_SNR.zip diff --git a/enpt/io/reader.py b/enpt/io/reader.py index 385af79..e3a6ce2 100644 --- a/enpt/io/reader.py +++ b/enpt/io/reader.py @@ -16,8 +16,7 @@ class L1B_Reader(object): self.logger = logger or logging.getLogger(__name__) self.cfg = user_inputs - @staticmethod - def read_inputdata(root_dir, observation_time: datetime, lon_lat_smpl=(15, 15), nsmile_coef=5, + def read_inputdata(self, root_dir, observation_time: datetime, lon_lat_smpl=(15, 15), nsmile_coef=5, snr_vnir=None, snr_swir=None): # TODO move to init? """Read L1B, DEM and spatial reference data. @@ -51,10 +50,12 @@ class L1B_Reader(object): # compute radiance and calculate snr L1_obj.DN2TOARadiance() + if snr_vnir is not None: + self.logger.info("Compute SNR for vnir: %s" % snr_vnir) + L1_obj.vnir.detector_meta.calc_snr_vnir(detector=L1_obj.vnir, snr_data_fn=snr_vnir) if snr_swir is not None: - L1_obj.vnir.detector_meta.calc_snr_vnir(data=L1_obj.vnir, snr_data_fn=snr_vnir) - if snr_swir is not None: - L1_obj.swir.detector_meta.calc_snr_swir(data=L1_obj.swir, snr_data_fn=snr_swir) + self.logger.info("Compute SNR for swir: %s" % snr_swir) + L1_obj.swir.detector_meta.calc_snr_swir(detector=L1_obj.swir, snr_data_fn=snr_swir) return L1_obj diff --git a/enpt/model/metadata.py b/enpt/model/metadata.py index 3d4dd33..60ae405 100644 --- a/enpt/model/metadata.py +++ b/enpt/model/metadata.py @@ -8,8 +8,6 @@ import numpy as np from scipy.interpolate import interp2d import spectral as sp -from ..model.images import EnMAP_Detector_SensorGeo - L1B_product_props = dict( xml_detector_label=dict( VNIR='detector1', @@ -124,28 +122,39 @@ class EnMAP_Metadata_L1B_Detector_SensorGeo(object): self.smile_coef # shape = (nwvl, nsmile_coef) ) # shape = (ncols, nwvl) - def calc_snr_vnir(self, data: EnMAP_Detector_SensorGeo, snr_data_fn: str): - """Compute EnMAP SNR from radiance data. + def calc_snr_vnir(self, detector, snr_data_fn: str) -> None: + """Compute EnMAP-VNIR SNR from radiance data. - :param data: Numpy array with radiance for scene + :param snr_data_fn: filename to SNR file + :param detector: vnir instance of EnMAP_Detector_SensorGeo """ - self.logger.info("Compute snr for: %s" % self.detector_name) + assert self.unit == 'mW m^-2 sr^-1 nm^-1' + self.logger.info("Compute snr for: %s using: %s" % (self.detector_name, snr_data_fn)) ds = sp.open_image(snr_data_fn) p_hg = ds[0:3, :, :] p_lg = ds[0:3, :, :] - l_th = ds[6, :, :] + l_th = np.squeeze(ds[6, :, :]) + self.snr = np.zeros(detector.data.shape) + for irow in range(self.nrows): + hg = detector.data[irow, :, :] < l_th + l_hg = detector.data[irow, hg] + self.snr[irow, :, :][hg] = p_hg[0, hg] + p_hg[1, hg] * l_hg + p_hg[2, hg] * l_hg**2 - return 500 * np.ones(data.shape, dtype=np.float) + lg = detector.data[irow, :, :] >= l_th + l_lg = detector.data[irow, lg] + self.snr[irow, :, :][lg] = p_lg[0, lg] + p_lg[1, lg] * l_lg + p_lg[2, lg] * l_lg ** 2 - def calc_snr_swir(self, data): + def calc_snr_swir(self, detector, snr_data_fn: str) -> None: """Compute EnMAP SNR from radiance data. - :param data: Numpy array with radiance for scene + :param snr_data_fn: filename to SNR file + :param detector: swir instance of EnMAP_Detector_SensorGeo """ - self.logger.info("Compute snr for: %s" % self.detector_name) - self.logger.warning("SNR model missing -> const. value of 500 is returned") - return 500 * np.ones(data.shape, dtype=np.float) + assert self.unit == 'mW m^-2 sr^-1 nm^-1' + self.logger.info("Compute snr for: %s using: %s" % (self.detector_name, snr_data_fn)) + pp = sp.open_image(snr_data_fn)[:, :, :] + self.snr = pp[0, :, :] + pp[1, :, :] * detector.data[:, :, :] + pp[2, :, :] * detector.data[:, :, :]**2 @staticmethod def interpolate_corners(ul: float, ur: float, ll: float, lr: float, nx: int, ny: int): diff --git a/tests/data/EnMAP_Sensor/EnMAP_Level_1B_SNR.zip b/tests/data/EnMAP_Sensor/EnMAP_Level_1B_SNR.zip new file mode 100644 index 0000000..a482a39 --- /dev/null +++ b/tests/data/EnMAP_Sensor/EnMAP_Level_1B_SNR.zip @@ -0,0 +1,3 @@ +version https://git-lfs.github.com/spec/v1 +oid sha256:4f32268975fd276091c41d8b256ad3b65653000cf99f0d94f34e1e91c912ea92 +size 13778 diff --git a/tests/test_l1b_reader.py b/tests/test_l1b_reader.py index b9585d8..4a45e63 100644 --- a/tests/test_l1b_reader.py +++ b/tests/test_l1b_reader.py @@ -21,6 +21,8 @@ class Test_L1B_Reader(unittest.TestCase): def setUp(self): self.pathList_testimages = glob(os.path.join(os.path.dirname(__file__), "data", "EnMAP_Level_1B", "*.zip")) + self.l1b_snr_file = glob(os.path.join(os.path.dirname(__file__), + "data", "EnMAP_Sensor", "EnMAP_Level_1B_SNR.zip"))[0] self.user_config = dict() def tearDown(self): @@ -34,14 +36,23 @@ class Test_L1B_Reader(unittest.TestCase): for l1b_file in self.pathList_testimages: with tempfile.TemporaryDirectory() as tmpdir: print("Tmp dir: %s" % tmpdir) - with zipfile.ZipFile(l1b_file, "r") as zf: + with zipfile.ZipFile(self.l1b_snr_file, "r") as zf: zf.extractall(tmpdir) - root_dir = os.path.join(tmpdir, os.listdir(tmpdir)[0]) - L1_obj = L1B_Reader(**self.user_config)\ - .read_inputdata(root_dir, observation_time=datetime(2015, 12, 7, 10)) + with zipfile.ZipFile(l1b_file, "r") as zf: + zf.extractall(tmpdir) - self.assertIsInstance(L1_obj, EnMAPL1Product_SensorGeo) + root_dir = os.path.join(tmpdir, os.path.basename(l1b_file).split(".zip")[0]) + L1_obj_no_snr = L1B_Reader(**self.user_config)\ + .read_inputdata(root_dir, observation_time=datetime(2015, 12, 7, 10)) + + L1_obj_with_snr = L1B_Reader(**self.user_config) \ + .read_inputdata(root_dir, observation_time=datetime(2015, 12, 7, 10), + snr_vnir=os.path.join(tmpdir, "SNR_D1.hdr"), + snr_swir=os.path.join(tmpdir, "SNR_D2.hdr")) + + self.assertIsInstance(L1_obj_no_snr, EnMAPL1Product_SensorGeo) + self.assertIsInstance(L1_obj_with_snr, EnMAPL1Product_SensorGeo) if __name__ == "__main__": -- GitLab From d3d85514b32477c6b3217381cb7684c69d5349d9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Hollstein?= Date: Thu, 31 Aug 2017 08:22:23 +0200 Subject: [PATCH 3/3] Add spectral library to requirements. --- setup.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/setup.py b/setup.py index 414d0d0..409dec2 100644 --- a/setup.py +++ b/setup.py @@ -9,8 +9,8 @@ with open('README.rst') as readme_file: with open('HISTORY.rst') as history_file: history = history_file.read() -requirements = [ - 'numpy', 'scipy', 'geoarray' # put package requirements here +requirements = [ # put package requirements here + 'numpy', 'scipy', 'geoarray', 'spectral' ] test_requirements = ['coverage', 'nose', 'nose-htmloutput', 'rednose'] -- GitLab