Skip to content

Update hypso.py #1

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
CameronLP wants to merge 4 commits into
base: master
Choose a base branch
from
Open

Update hypso.py #1

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
a8ed42d
Update hypso.py
CameronLP Feb 5, 2026
d012426
Update Level1_HYPSO class
CameronLP Feb 13, 2026
d53df19
Fix SyntaxError: f-string: unmatched '[' error
CameronLP Feb 13, 2026
e6ae1db
Fix line 26 SyntaxError
CameronLP Feb 13, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions eoread/autodetect.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@ def Level1(path: Path, **kwargs) -> xr.Dataset:
'Path does not correspond to a level 1 product'

# Import reader and read provided file
assert dict_pattern['reader'] != '', f'No reader exists for {dict_pattern['Name']}'
assert dict_pattern['reader'] != '', f"No reader exists for {dict_pattern['Name']}"
to_import = dict_pattern['reader'].split()
module = importlib.import_module(to_import[0])
try: reader = getattr(module, to_import[1].replace('Level','Level1'))
Expand All @@ -42,7 +42,7 @@ def Level2(path: Path, **kwargs) -> xr.Dataset:
'Path does not correspond to a level 2 product'

# Import reader and read provided file
assert dict_pattern['reader'] != '', f'No reader exists for {dict_pattern['Name']}'
assert dict_pattern['reader'] != '', f"No reader exists for {dict_pattern['Name']}"
to_import = dict_pattern['reader'].split()
module = importlib.import_module(to_import[0])
try: reader = getattr(module, to_import[1].replace('Level','Level2'))
Expand Down
140 changes: 119 additions & 21 deletions eoread/hypso.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3,16 +3,24 @@
import xarray as xr

from eoread.utils import filter_metadata
# from eotools.solar_irradiance import solar_irradiance
from eotools.solar_irradiance import solar_irradiance
from core.interpolate import interp, Linear
from core.geo import n
from core.tools import drop_unused_dims
from core import env, log



#rtoa = Var("Rtoa", attrs={'desc': "Top of Atmosphere reflectance"}, dtype='float32')





def Level1_HYPSO(filepath: str|Path,
chunks: int|tuple = 500,
metadata_template: list = None,
v1_compat: bool = False) -> xr.Dataset:
v1_compat: bool = True) -> xr.Dataset:
'''
Read an NTNU HYPSO Level1 product as an xarray.Dataset
Formats the Dataset so that it contains the TOA radiances,
Expand All @@ -29,10 +37,37 @@ def Level1_HYPSO(filepath: str|Path,
filepath = Path(filepath)
assert filepath.exists(), 'File does not exists'


ds_root = xr.open_datatree(filepath, engine='h5netcdf')

match str(ds_root.attrs['processing_level']).upper():

case "L1C":
hypso_product_name = "Lt"
polymer_product_name = n.ltoa
units = 'W/sr/m^2'
load_f0 = True

case "L1D":
hypso_product_name = "Rhot"
polymer_product_name = n.rtoa
units = ''
load_f0 = False

case _:
log.debug('Unsupported HYPSO product level.')
return


print(ds_root['metadata']['corrections'].attrs['radiometric_coefficents_version'])


if isinstance(chunks, int): chunks = [chunks]*2
ds_products = ds_root["products"].to_dataset()
ds_nav = ds_root["navigation"].to_dataset()
try:
ds_nav = ds_root["navigation"].to_dataset()
except KeyError:
ds_nav = ds_root["geometry"].to_dataset()

# get _indirect geographical coordinates and angles if available
log.debug('Read and compute geometric angles')
Expand All @@ -44,23 +79,72 @@ def Level1_HYPSO(filepath: str|Path,
ds[str(n.saa)] = ds_nav["solar_azimuth"].chunk(chunks)

log.debug('Read top of atmosphere data')
ds[str(n.ltoa)] = ds_products['Lt'].chunk(list(chunks)+[1])
ds = ds.rename(lines=str(n.rows), samples=str(n.columns), bands=str(n.bands))
ds[str(n.ltoa)].attrs['unit'] = 'W/sr/m^2'

try:
ds[str(polymer_product_name)] = ds_products[hypso_product_name].chunk(list(chunks)+[1])
ds = ds.rename(lines=str(n.rows), samples=str(n.columns), bands=str(n.bands))
ds[str(polymer_product_name)].attrs['unit'] = 'W/sr/m^2'
wavelengths = ds_products[hypso_product_name].wavelengths
except KeyError:

import numpy as np

wavelengths = []

for band_idx, band_key in enumerate(ds_products.keys()):

if band_idx == 0:
rows, columns = ds_products[band_key].shape
bands = len(ds_products.keys())
datacube = np.empty((rows,columns,bands))

band_data = np.array(ds_products[band_key][:], dtype='float')

wave = ds_products[band_key].attrs["wavelength"]

wavelengths.append(wave)

datacube[:,:,band_idx] = band_data


ds[str(polymer_product_name)] = xr.DataArray(datacube, name=hypso_product_name, dims=['lines', 'samples', 'bands']).chunk(list(chunks)+[1])
ds = ds.rename(lines=str(n.rows), samples=str(n.columns), bands=str(n.bands))
ds[str(polymer_product_name)].attrs['unit'] = units

wavelengths = np.array(wavelengths)

wavelengths = np.around(wavelengths,1) # round to one decimal like hypso-package L1 processing
wavelengths = np.array(wavelengths).astype(np.int32) # convert to int like hypso-package L1 processing
wavelengths = np.array(wavelengths, dtype='float32') # convert to float for NetCDF writing compatibility



log.debug('Extract central wavelength')
ds = ds.assign_coords({str(n.bands): da.arange(len(ds[str(n.bands)]))+1})
ds = ds.assign({str(n.cwav): ((str(n.bands)), ds_products['Lt'].wavelengths),
#ds = ds.assign_coords({str(n.bands): da.arange(len(ds[str(n.bands)]))+1})
ds = ds.assign_coords({str(n.bands): wavelengths})

ds = ds.assign({str(n.cwav): ((str(n.bands)), wavelengths),
str(n.wav): ((str(n.bands)), wavelengths),
str(n.bnames): ((str(n.bands)), ds[str(n.bands)].data.astype(str))})

# # read solar irradiance
# F0 = solar_irradiance("LISIRD", variant="1nm")
# ds["F0"] = interp(F0, wavelength=Linear(ds.wav))
# # convert it to a unit compatible with Ltoa
# assert ds.F0.units == "W m-2 nm-1"
# ds["F0"] = ds.F0 * 1000
# ds.F0.attrs.update(units="W m-2 um-1")
# ds = ds.rename(lines="y", samples="x")
if load_f0:
# read solar irradiance
F0 = solar_irradiance("LISIRD", variant="1nm")

print(F0)
interpolated_F0 = interp(F0, wavelength=Linear(wavelengths))
#ds["F0"] = interp(F0, wavelength=Linear(ds.wav))

# convert it to a unit compatible with Ltoa
#assert ds.F0.units == "W m-2 nm-1"
#ds["F0"] = ds.F0 * 1000
#ds.F0.attrs.update(units="W m-2 um-1")
#ds = ds.rename(lines="y", samples="x")


interpolated_F0 = interpolated_F0*1000
ds["F0"] = xr.DataArray(interpolated_F0, name="F0")
#ds.rename(lines=str(n.rows), samples=str(n.columns))
ds["F0"].attrs.update(units="W m-2 um-1")

# acquisition datetime
log.debug('Add important attributes')
Expand All @@ -69,11 +153,25 @@ def Level1_HYPSO(filepath: str|Path,
ds.attrs[str(n.resolution)] = 40
ds.attrs[str(n.product_name)] = filepath.name
ds.attrs[str(n.input_directory)] = str(filepath.parent)
ds.attrs[str(n.datetime)] = ds_root.attrs['date_aquired']


try:
ds.attrs[str(n.datetime)] = ds_root.attrs['date_aquired']
except KeyError:

timestr = ds_root.attrs['timestamp_acquired']
if timestr.endswith('+00:00Z'):
timestr = timestr.replace('+00:00Z', 'Z') # or '+00:00'


ds.attrs[str(n.datetime)] = timestr



filter_fn = (lambda x,y: x) if metadata_template is None else filter_metadata
ds.attrs['metadata'] = filter_fn(ds_root.attrs, metadata_template)

ds.attrs.pop("metadata")

# ds[naming.flags] = xr.zeros_like(ds.vza, dtype=naming.flags_dtype)

if v1_compat: return _v1_compat(ds)
Expand All @@ -95,6 +193,6 @@ def get_sample(level: int=1) -> Path:
def _v1_compat(ds):

# Add flags
ds["flags"] = xr.zeros_like(ds.vza, dtype="uint8")
ds["flags"] = xr.zeros_like(ds.vza, dtype="uint16")

return ds
return ds