Commit ea78f442 authored by Daniel Scheffler's avatar Daniel Scheffler
Browse files

Revised SRF object.

Former-commit-id: df54bc8b
Former-commit-id: e7188b72
parent ba535300
......@@ -56,7 +56,7 @@ class L2B_object(L2A_object):
orig_CWLs, target_CWLs = np.array(source_CWLs), np.array(target_CWLs)
outarr = interp1d(np.array(orig_CWLs), arrcube, axis=2, kind='linear', fill_value='extrapolate')(target_CWLs)
outarr = outarr.astype(np.int16)
assert outarr.shape == (arrcube.shape[:2] + (len(target_CWLs), ))
assert outarr.shape == tuple([*arrcube.shape[:2], len(target_CWLs)])
return outarr
......
......@@ -10,7 +10,9 @@ import tarfile
import warnings
import zipfile
from tempfile import NamedTemporaryFile as tempFile
from typing import TYPE_CHECKING
from logging import Logger
from matplotlib import pyplot as plt
from typing import Dict # noqa F401 # flake8 issue
import dill
import numpy as np
......@@ -32,9 +34,6 @@ from geoarray import GeoArray
from py_tools_ds.geo.coord_calc import corner_coord_to_minmax
from py_tools_ds.geo.coord_trafo import transform_any_prj
if TYPE_CHECKING:
from logging import Logger # noqa F401 # flake8 issue
def read_ENVIfile(path, arr_shape, arr_pos, logger=None, return_meta=False, q=0):
hdr_path = os.path.splitext(path)[0] + '.hdr' if not os.path.splitext(path)[1] == '.hdr' else path
......@@ -184,7 +183,8 @@ def SRF_reader(GMS_identifier, v=False):
:param v: verbose mode
"""
satellite, sensor, logger = (GMS_identifier['Satellite'], GMS_identifier['Sensor'], GMS_identifier['logger'])
satellite, sensor = GMS_identifier['Satellite'], GMS_identifier['Sensor']
logger = GMS_identifier['logger'] or Logger(__name__)
LayerBandsAssignment = META.get_LayerBandsAssignment(GMS_identifier)
SRF_dict = collections.OrderedDict()
SRF_dir = PG.get_path_srf_file(GMS_identifier)
......@@ -221,38 +221,20 @@ class SRF(object):
raise ValueError('Unknown wavelength unit %s.' % wvl_unit)
self.srfs_wvl = []
self.srfs = []
self.srfs = {}
self.srfs_norm01 = {}
self.bands = []
self.wvl = []
self.wvl_unit = wvl_unit
self.format_bandnames = format_bandnames
self.satellite = GMS_identifier['Satellite'] if GMS_identifier else ''
self.sensor = GMS_identifier['Sensor'] if GMS_identifier else ''
self.conv = {}
self.v = v
if GMS_identifier:
self.from_GMS_identifier(GMS_identifier)
def _set_attrs(self, wvl, srfs):
for bN in srfs:
srf = np.array(srfs[bN], dtype=np.float)
srfs[bN] = srf / np.trapz(x=wvl, y=srf) # TODO seems like we NEED nanometers here; BUT WHY??
self.srfs_wvl = np.array(wvl)
self.srfs = srfs
self.bands = sorted(list(self.srfs.keys()))
# FIXME this is not the GMS algorithm to calculate center wavelengths
# calculate center wavelengths
# TODO seems like we NEED nanometers here; BUT WHY??:
self.wvl = np.array([np.trapz(x=self.srfs_wvl, y=self.srfs_wvl * self.srfs[band]) for band in self.bands])
# self.wvl = self.wvl if self.wvl_unit == 'micrometers' else np.array([int(i) for i in self.wvl])
self.srfs_wvl = self.srfs_wvl if self.wvl_unit == 'nanometers' else self.srfs_wvl / 1000
self.wvl = self.wvl if self.wvl_unit == 'nanometers' else self.wvl / 1000
self.conv = {}
self.conv.update({key: value for key, value in zip(self.bands, self.wvl)})
self.conv.update({value: key for key, value in zip(self.bands, self.wvl)})
def from_GMS_identifier(self, GMS_identifier):
srf_dict = SRF_reader(GMS_identifier, v=self.v)
return self.from_dict(srf_dict)
......@@ -278,7 +260,30 @@ class SRF(object):
df = df.join(Series(srfs, index=wvls, name=bandname))
df = df.fillna(0)
wvl = np.array(df.index.astype(np.float))
self._set_attrs(wvl=wvl, srfs=df.to_dict(orient='series'))
srfs_norm01 = df.to_dict(orient='series') # type: Dict[Series]
################
# set attributes
################
for bN in srfs_norm01:
self.srfs_norm01[bN] = srf = np.array(srfs_norm01[bN], dtype=np.float)
self.srfs[bN] = srf / np.trapz(x=wvl, y=srf) # TODO seems like we NEED nanometers here; BUT WHY??
self.srfs_wvl = np.array(wvl)
self.bands = sorted(list(self.srfs.keys()))
# FIXME this is not the GMS algorithm to calculate center wavelengths
# calculate center wavelengths
# TODO seems like we NEED nanometers here; BUT WHY??:
self.wvl = np.array([np.trapz(x=self.srfs_wvl, y=self.srfs_wvl * self.srfs[band]) for band in self.bands])
# self.wvl = self.wvl if self.wvl_unit == 'micrometers' else np.array([int(i) for i in self.wvl])
self.srfs_wvl = self.srfs_wvl if self.wvl_unit == 'nanometers' else self.srfs_wvl / 1000
self.wvl = self.wvl if self.wvl_unit == 'nanometers' else self.wvl / 1000
self.conv.update({key: value for key, value in zip(self.bands, self.wvl)})
self.conv.update({value: key for key, value in zip(self.bands, self.wvl)})
return self
......@@ -303,6 +308,22 @@ class SRF(object):
def __getitem__(self, band):
return self.srfs[band]
def plot_srfs(self, figsize=(15, 5)):
# type: (tuple) -> None
"""Show a plot of all spectral response functions.
:param figsize: figue size of the plot
"""
from ..misc.helper_functions import sorted_nicely
plt.figure(figsize=figsize)
for band in sorted_nicely(self.bands):
plt.plot(self.srfs_wvl, list(self.srfs_norm01[band]), '-', label='Band %s' % band)
plt.title(' '.join(['Spectral response functions', self.satellite, self.sensor]))
plt.xlabel('wavelength [%s]' % self.wvl_unit)
plt.ylabel('spectral response')
plt.legend(loc='upper right')
def pickle_SRF_DB(L1A_Instances):
list_GMS_identifiers = [i.GMS_identifier for i in L1A_Instances]
......
py_tools_ds>=0.9.1
geoarray>=0.6.12
arosics>0.6.2
arosics>=0.6.2
git+https://gitext.gfz-potsdam.de/EnMAP/sicor.git
matplotlib
numpy
......
Markdown is supported
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