From cef002a57fc0a96afc95e5bb98276c89877dfde8 Mon Sep 17 00:00:00 2001 From: David Hoese Date: Tue, 30 Nov 2021 20:40:32 -0600 Subject: [PATCH 1/6] Refactor SZA and cos(SZA) generation to reduce duplicate computations --- satpy/modifiers/angles.py | 105 +++++++++++++++++++++++++++++++----- satpy/modifiers/geometry.py | 32 ++++++----- satpy/utils.py | 36 ------------- 3 files changed, 113 insertions(+), 60 deletions(-) diff --git a/satpy/modifiers/angles.py b/satpy/modifiers/angles.py index fb06431c9a..c51d34ea94 100644 --- a/satpy/modifiers/angles.py +++ b/satpy/modifiers/angles.py @@ -27,10 +27,11 @@ from typing import Any, Callable, Optional, Union import dask -import dask.array as da import numpy as np import xarray as xr -from pyorbital.astronomy import get_alt_az, sun_zenith_angle +from dask import array as da +from pyorbital.astronomy import cos_zen as pyob_cos_zen +from pyorbital.astronomy import get_alt_az from pyorbital.orbital import get_observer_look from pyresample.geometry import AreaDefinition, SwathDefinition @@ -253,7 +254,7 @@ def get_angles(data_arr: xr.DataArray) -> tuple[xr.DataArray, xr.DataArray, xr.D Returns: Four DataArrays representing sensor azimuth angle, sensor zenith angle, - solar azimuth angle, and solar zenith angle. + solar azimuth angle, and solar zenith angle. All values are in degrees. """ sata, satz = _get_sensor_angles(data_arr) @@ -274,6 +275,18 @@ def get_satellite_zenith_angle(data_arr: xr.DataArray) -> xr.DataArray: return satz +def get_cos_sza(data_arr: xr.DataArray) -> xr.DataArray: + """Generate the cosine of the solar zenith angle for the provided data. + + Returns: + DataArray with the same shape as ``data_arr``. Units are radians. + + """ + lons, lats = _get_valid_lonlats(data_arr.attrs["area"], data_arr.chunks) + cos_sza = _get_cos_sza(data_arr.attrs["start_time"], lons, lats) + return _geo_dask_to_data_array(cos_sza) + + @cache_to_zarr_if("cache_lonlats") def _get_valid_lonlats(area: PRGeometry, chunks: Union[int, str] = "auto") -> tuple[da.Array, da.Array]: with ignore_invalid_float_warnings(): @@ -285,21 +298,43 @@ def _get_valid_lonlats(area: PRGeometry, chunks: Union[int, str] = "auto") -> tu def _get_sun_angles(data_arr: xr.DataArray) -> tuple[xr.DataArray, xr.DataArray]: lons, lats = _get_valid_lonlats(data_arr.attrs["area"], data_arr.data.chunks) - res = da.map_blocks(_get_sun_angles_wrapper, lons, lats, - data_arr.attrs["start_time"], - dtype=lons.dtype, meta=np.array((), dtype=lons.dtype), - new_axis=[0], chunks=(2,) + lons.chunks) - suna = _geo_dask_to_data_array(res[0]) - sunz = _geo_dask_to_data_array(res[1]) + suna = da.map_blocks(_get_alt_az_wrapper, lons, lats, + data_arr.attrs["start_time"], + dtype=lons.dtype, meta=np.array((), dtype=lons.dtype), + chunks=lons.chunks) + cos_sza = _get_cos_sza(data_arr.attrs["start_time"], lons, lats) + sunz = np.rad2deg(np.arccos(cos_sza)) + suna = _geo_dask_to_data_array(suna) + sunz = _geo_dask_to_data_array(sunz) return suna, sunz -def _get_sun_angles_wrapper(lons: da.Array, lats: da.Array, start_time: datetime) -> np.ndarray: +def _get_cos_sza(utc_time, lons, lats): + cos_sza = da.map_blocks(_cos_zen_wrapper, + lons, lats, utc_time, + meta=np.array((), dtype=lons.dtype), + dtype=lons.dtype, + chunks=lons.chunks) + return cos_sza + + +def _cos_zen_wrapper(lons, lats, utc_time): + with ignore_invalid_float_warnings(): + return pyob_cos_zen(utc_time, lons, lats) + + +def _get_alt_az_wrapper(lons: np.ndarray, lats: np.ndarray, start_time: datetime) -> np.ndarray: with ignore_invalid_float_warnings(): suna = get_alt_az(start_time, lons, lats)[1] suna = np.rad2deg(suna) - sunz = sun_zenith_angle(start_time, lons, lats) - return np.stack([suna, sunz]) + return suna + + +# def _get_sza_wrapper(lons: np.ndarray, lats: np.ndarray, start_time: datetime) -> np.ndarray: +# with ignore_invalid_float_warnings(): +# # return np.rad2deg(np.arccos(cos_zen(utc_time, lon, lat))) +# sunz = sun_zenith_angle(start_time, lons, lats) +# return sunz def _get_sensor_angles(data_arr: xr.DataArray) -> tuple[xr.DataArray, xr.DataArray]: @@ -332,3 +367,49 @@ def _get_sensor_angles_wrapper(lons, lats, start_time, sat_lon, sat_lat, sat_alt lons, lats, 0) satz = 90 - satel return np.stack([sata, satz]) + + +def sunzen_corr_cos(data, cos_zen, limit=88., max_sza=95.): + """Perform Sun zenith angle correction. + + The correction is based on the provided cosine of the zenith + angle (``cos_zen``). The correction is limited + to ``limit`` degrees (default: 88.0 degrees). For larger zenith + angles, the correction is the same as at the ``limit`` if ``max_sza`` + is `None`. The default behavior is to gradually reduce the correction + past ``limit`` degrees up to ``max_sza`` where the correction becomes + 0. Both ``data`` and ``cos_zen`` should be 2D arrays of the same shape. + + """ + return da.map_blocks(_sunzen_corr_cos_wrapper, + data, cos_zen, limit, max_sza, + meta=np.array((), dtype=data.dtype), + chunks=data.chunks) + + +def _sunzen_corr_cos_wrapper(data, cos_zen, limit, max_sza): + # Convert the zenith angle limit to cosine of zenith angle + limit_rad = np.deg2rad(limit) + limit_cos = np.cos(limit_rad) + max_sza_rad = np.deg2rad(max_sza) if max_sza is not None else max_sza + + # Cosine correction + corr = 1. / cos_zen + if max_sza is not None: + # gradually fall off for larger zenith angle + grad_factor = (np.arccos(cos_zen) - limit_rad) / (max_sza_rad - limit_rad) + # invert the factor so maximum correction is done at `limit` and falls off later + grad_factor = 1. - np.log(grad_factor + 1) / np.log(2) + # make sure we don't make anything negative + grad_factor = grad_factor.clip(0.) + else: + # Use constant value (the limit) for larger zenith angles + grad_factor = 1. + # corr[cos_zen <= limit_cos] = grad_factor / limit_cos + # corr = corr.where(cos_zen > limit_cos, grad_factor / limit_cos) + corr = np.where(cos_zen > limit_cos, corr, grad_factor / limit_cos) + # Force "night" pixels to 0 (where SZA is invalid) + corr[np.isnan(cos_zen)] = 0 + # corr = corr.where(cos_zen.notnull(), 0) + + return data * corr diff --git a/satpy/modifiers/geometry.py b/satpy/modifiers/geometry.py index 0913fb42d1..6188ca7b05 100644 --- a/satpy/modifiers/geometry.py +++ b/satpy/modifiers/geometry.py @@ -27,7 +27,8 @@ import xarray as xr from satpy.modifiers import ModifierBase -from satpy.utils import atmospheric_path_length_correction, sunzen_corr_cos +from satpy.modifiers.angles import sunzen_corr_cos +from satpy.utils import atmospheric_path_length_correction logger = logging.getLogger(__name__) @@ -63,17 +64,21 @@ def __call__(self, projectables, **info): logger.debug("Applying sun zen correction") coszen = self.coszen_cache.get(key) if coszen is None and not info.get('optional_datasets'): - # we were not given SZA, generate SZA then calculate cos(SZA) - from pyorbital.astronomy import cos_zen + # we were not given SZA, generate cos(SZA) logger.debug("Computing sun zenith angles.") - lons, lats = vis.attrs["area"].get_lonlats(chunks=vis.data.chunks) - - coords = {} - if 'y' in vis.coords and 'x' in vis.coords: - coords['y'] = vis['y'] - coords['x'] = vis['x'] - coszen = xr.DataArray(cos_zen(vis.attrs["start_time"], lons, lats), - dims=['y', 'x'], coords=coords) + # coords = {} + # if 'y' in vis.coords and 'x' in vis.coords: + # coords['y'] = vis['y'] + # coords['x'] = vis['x'] + + from .angles import get_cos_sza + coszen = get_cos_sza(vis) + # lons, lats = vis.attrs["area"].get_lonlats(chunks=vis.data.chunks) + # coszen = xr.DataArray(cos_zen(vis.attrs["start_time"], lons, lats), + # dims=['y', 'x'], coords=coords) + # import dask.array as da + # coszen = xr.DataArray(da.ones_like(lons), + # dims=['y', 'x'], coords=coords) if self.max_sza is not None: coszen = coszen.where(coszen >= self.max_sza_cos) self.coszen_cache[key] = coszen @@ -130,7 +135,10 @@ def __init__(self, correction_limit=88., **kwargs): def _apply_correction(self, proj, coszen): logger.debug("Apply the standard sun-zenith correction [1/cos(sunz)]") - return sunzen_corr_cos(proj, coszen, limit=self.correction_limit, max_sza=self.max_sza) + print("Applying sunzen: ", proj.chunks == coszen.chunks) + res = proj.copy() + res.data = sunzen_corr_cos(proj.data, coszen.data, limit=self.correction_limit, max_sza=self.max_sza) + return res class EffectiveSolarPathLengthCorrector(SunZenithCorrectorBase): diff --git a/satpy/utils.py b/satpy/utils.py index 7c362a672d..af61b4ff51 100644 --- a/satpy/utils.py +++ b/satpy/utils.py @@ -240,42 +240,6 @@ def _get_sunz_corr_li_and_shibata(cos_zen): return 24.35 / (2. * cos_zen + np.sqrt(498.5225 * cos_zen**2 + 1)) -def sunzen_corr_cos(data, cos_zen, limit=88., max_sza=95.): - """Perform Sun zenith angle correction. - - The correction is based on the provided cosine of the zenith - angle (``cos_zen``). The correction is limited - to ``limit`` degrees (default: 88.0 degrees). For larger zenith - angles, the correction is the same as at the ``limit`` if ``max_sza`` - is `None`. The default behavior is to gradually reduce the correction - past ``limit`` degrees up to ``max_sza`` where the correction becomes - 0. Both ``data`` and ``cos_zen`` should be 2D arrays of the same shape. - - """ - # Convert the zenith angle limit to cosine of zenith angle - limit_rad = np.deg2rad(limit) - limit_cos = np.cos(limit_rad) - max_sza_rad = np.deg2rad(max_sza) if max_sza is not None else max_sza - - # Cosine correction - corr = 1. / cos_zen - if max_sza is not None: - # gradually fall off for larger zenith angle - grad_factor = (np.arccos(cos_zen) - limit_rad) / (max_sza_rad - limit_rad) - # invert the factor so maximum correction is done at `limit` and falls off later - grad_factor = 1. - np.log(grad_factor + 1) / np.log(2) - # make sure we don't make anything negative - grad_factor = grad_factor.clip(0.) - else: - # Use constant value (the limit) for larger zenith angles - grad_factor = 1. - corr = corr.where(cos_zen > limit_cos, grad_factor / limit_cos) - # Force "night" pixels to 0 (where SZA is invalid) - corr = corr.where(cos_zen.notnull(), 0) - - return data * corr - - def atmospheric_path_length_correction(data, cos_zen, limit=88., max_sza=95.): """Perform Sun zenith angle correction. From 85dfa745b6c2895d968409780c650e783328899f Mon Sep 17 00:00:00 2001 From: David Hoese Date: Tue, 30 Nov 2021 20:51:46 -0600 Subject: [PATCH 2/6] Refactor SZA and cos(SZA) generation to reduce duplicate computations --- satpy/modifiers/geometry.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/satpy/modifiers/geometry.py b/satpy/modifiers/geometry.py index 6188ca7b05..ab2df7de88 100644 --- a/satpy/modifiers/geometry.py +++ b/satpy/modifiers/geometry.py @@ -17,6 +17,8 @@ # satpy. If not, see . """Modifier classes for corrections based on sun and other angles.""" +from __future__ import annotations + import logging import time from datetime import datetime @@ -135,7 +137,6 @@ def __init__(self, correction_limit=88., **kwargs): def _apply_correction(self, proj, coszen): logger.debug("Apply the standard sun-zenith correction [1/cos(sunz)]") - print("Applying sunzen: ", proj.chunks == coszen.chunks) res = proj.copy() res.data = sunzen_corr_cos(proj.data, coszen.data, limit=self.correction_limit, max_sza=self.max_sza) return res From 19c916b4e762f6061d4cc2ae578ec55b5039a688 Mon Sep 17 00:00:00 2001 From: David Hoese Date: Wed, 1 Dec 2021 09:43:48 -0600 Subject: [PATCH 3/6] Remove commented code and other code cleanup --- satpy/modifiers/angles.py | 35 +++++++++++++++-------------------- satpy/modifiers/geometry.py | 11 ----------- 2 files changed, 15 insertions(+), 31 deletions(-) diff --git a/satpy/modifiers/angles.py b/satpy/modifiers/angles.py index c51d34ea94..0d30875e5c 100644 --- a/satpy/modifiers/angles.py +++ b/satpy/modifiers/angles.py @@ -298,7 +298,7 @@ def _get_valid_lonlats(area: PRGeometry, chunks: Union[int, str] = "auto") -> tu def _get_sun_angles(data_arr: xr.DataArray) -> tuple[xr.DataArray, xr.DataArray]: lons, lats = _get_valid_lonlats(data_arr.attrs["area"], data_arr.data.chunks) - suna = da.map_blocks(_get_alt_az_wrapper, lons, lats, + suna = da.map_blocks(_get_sun_azimuth_ndarray, lons, lats, data_arr.attrs["start_time"], dtype=lons.dtype, meta=np.array((), dtype=lons.dtype), chunks=lons.chunks) @@ -310,7 +310,7 @@ def _get_sun_angles(data_arr: xr.DataArray) -> tuple[xr.DataArray, xr.DataArray] def _get_cos_sza(utc_time, lons, lats): - cos_sza = da.map_blocks(_cos_zen_wrapper, + cos_sza = da.map_blocks(_cos_zen_ndarray, lons, lats, utc_time, meta=np.array((), dtype=lons.dtype), dtype=lons.dtype, @@ -318,25 +318,18 @@ def _get_cos_sza(utc_time, lons, lats): return cos_sza -def _cos_zen_wrapper(lons, lats, utc_time): +def _cos_zen_ndarray(lons, lats, utc_time): with ignore_invalid_float_warnings(): return pyob_cos_zen(utc_time, lons, lats) -def _get_alt_az_wrapper(lons: np.ndarray, lats: np.ndarray, start_time: datetime) -> np.ndarray: +def _get_sun_azimuth_ndarray(lons: np.ndarray, lats: np.ndarray, start_time: datetime) -> np.ndarray: with ignore_invalid_float_warnings(): suna = get_alt_az(start_time, lons, lats)[1] suna = np.rad2deg(suna) return suna -# def _get_sza_wrapper(lons: np.ndarray, lats: np.ndarray, start_time: datetime) -> np.ndarray: -# with ignore_invalid_float_warnings(): -# # return np.rad2deg(np.arccos(cos_zen(utc_time, lon, lat))) -# sunz = sun_zenith_angle(start_time, lons, lats) -# return sunz - - def _get_sensor_angles(data_arr: xr.DataArray) -> tuple[xr.DataArray, xr.DataArray]: sat_lon, sat_lat, sat_alt = get_satpos(data_arr) area_def = data_arr.attrs["area"] @@ -351,13 +344,13 @@ def _get_sensor_angles(data_arr: xr.DataArray) -> tuple[xr.DataArray, xr.DataArr @cache_to_zarr_if("cache_sensor_angles", sanitize_args_func=_sanitize_observer_look_args) def _get_sensor_angles_from_sat_pos(sat_lon, sat_lat, sat_alt, start_time, area_def, chunks): lons, lats = _get_valid_lonlats(area_def, chunks) - res = da.map_blocks(_get_sensor_angles_wrapper, lons, lats, start_time, sat_lon, sat_lat, sat_alt, + res = da.map_blocks(_get_sensor_angles_ndarray, lons, lats, start_time, sat_lon, sat_lat, sat_alt, dtype=lons.dtype, meta=np.array((), dtype=lons.dtype), new_axis=[0], chunks=(2,) + lons.chunks) return res[0], res[1] -def _get_sensor_angles_wrapper(lons, lats, start_time, sat_lon, sat_lat, sat_alt) -> np.ndarray: +def _get_sensor_angles_ndarray(lons, lats, start_time, sat_lon, sat_lat, sat_alt) -> np.ndarray: with ignore_invalid_float_warnings(): sata, satel = get_observer_look( sat_lon, @@ -369,7 +362,10 @@ def _get_sensor_angles_wrapper(lons, lats, start_time, sat_lon, sat_lat, sat_alt return np.stack([sata, satz]) -def sunzen_corr_cos(data, cos_zen, limit=88., max_sza=95.): +def sunzen_corr_cos(data: da.Array, + cos_zen: da.Array, + limit: float = 88., + max_sza: Optional[float] = 95.) -> da.Array: """Perform Sun zenith angle correction. The correction is based on the provided cosine of the zenith @@ -381,13 +377,16 @@ def sunzen_corr_cos(data, cos_zen, limit=88., max_sza=95.): 0. Both ``data`` and ``cos_zen`` should be 2D arrays of the same shape. """ - return da.map_blocks(_sunzen_corr_cos_wrapper, + return da.map_blocks(_sunzen_corr_cos_ndarray, data, cos_zen, limit, max_sza, meta=np.array((), dtype=data.dtype), chunks=data.chunks) -def _sunzen_corr_cos_wrapper(data, cos_zen, limit, max_sza): +def _sunzen_corr_cos_ndarray(data: np.ndarray, + cos_zen: np.ndarray, + limit: float, + max_sza: Optional[float]) -> np.ndarray: # Convert the zenith angle limit to cosine of zenith angle limit_rad = np.deg2rad(limit) limit_cos = np.cos(limit_rad) @@ -405,11 +404,7 @@ def _sunzen_corr_cos_wrapper(data, cos_zen, limit, max_sza): else: # Use constant value (the limit) for larger zenith angles grad_factor = 1. - # corr[cos_zen <= limit_cos] = grad_factor / limit_cos - # corr = corr.where(cos_zen > limit_cos, grad_factor / limit_cos) corr = np.where(cos_zen > limit_cos, corr, grad_factor / limit_cos) # Force "night" pixels to 0 (where SZA is invalid) corr[np.isnan(cos_zen)] = 0 - # corr = corr.where(cos_zen.notnull(), 0) - return data * corr diff --git a/satpy/modifiers/geometry.py b/satpy/modifiers/geometry.py index ab2df7de88..8bbd1fd47d 100644 --- a/satpy/modifiers/geometry.py +++ b/satpy/modifiers/geometry.py @@ -68,19 +68,8 @@ def __call__(self, projectables, **info): if coszen is None and not info.get('optional_datasets'): # we were not given SZA, generate cos(SZA) logger.debug("Computing sun zenith angles.") - # coords = {} - # if 'y' in vis.coords and 'x' in vis.coords: - # coords['y'] = vis['y'] - # coords['x'] = vis['x'] - from .angles import get_cos_sza coszen = get_cos_sza(vis) - # lons, lats = vis.attrs["area"].get_lonlats(chunks=vis.data.chunks) - # coszen = xr.DataArray(cos_zen(vis.attrs["start_time"], lons, lats), - # dims=['y', 'x'], coords=coords) - # import dask.array as da - # coszen = xr.DataArray(da.ones_like(lons), - # dims=['y', 'x'], coords=coords) if self.max_sza is not None: coszen = coszen.where(coszen >= self.max_sza_cos) self.coszen_cache[key] = coszen From 309e8f741fa0af211d195be9c78ac63508e2b84e Mon Sep 17 00:00:00 2001 From: David Hoese Date: Wed, 1 Dec 2021 13:48:06 -0600 Subject: [PATCH 4/6] Fix type annotation for chunks argument _get_valid_lonlats --- satpy/modifiers/angles.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/satpy/modifiers/angles.py b/satpy/modifiers/angles.py index 0d30875e5c..3f57e070a2 100644 --- a/satpy/modifiers/angles.py +++ b/satpy/modifiers/angles.py @@ -288,7 +288,7 @@ def get_cos_sza(data_arr: xr.DataArray) -> xr.DataArray: @cache_to_zarr_if("cache_lonlats") -def _get_valid_lonlats(area: PRGeometry, chunks: Union[int, str] = "auto") -> tuple[da.Array, da.Array]: +def _get_valid_lonlats(area: PRGeometry, chunks: Union[int, str, tuple] = "auto") -> tuple[da.Array, da.Array]: with ignore_invalid_float_warnings(): lons, lats = area.get_lonlats(chunks=chunks) lons = da.where(lons >= 1e30, np.nan, lons) From e969751e262e7040e0eb4ceaf1fb951b3145cd7e Mon Sep 17 00:00:00 2001 From: David Hoese Date: Thu, 2 Dec 2021 21:05:16 -0600 Subject: [PATCH 5/6] Add a few additional cleanups to dask usage --- satpy/modifiers/atmosphere.py | 2 -- satpy/resample.py | 6 ++++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/satpy/modifiers/atmosphere.py b/satpy/modifiers/atmosphere.py index 431b303f40..ad33505747 100644 --- a/satpy/modifiers/atmosphere.py +++ b/satpy/modifiers/atmosphere.py @@ -45,11 +45,9 @@ def __call__(self, projectables, optional_datasets=None, **info): if not optional_datasets or len(optional_datasets) != 4: vis, red = self.match_data_arrays(projectables) sata, satz, suna, sunz = get_angles(vis) - red.data = da.rechunk(red.data, vis.data.chunks) else: vis, red, sata, satz, suna, sunz = self.match_data_arrays( projectables + optional_datasets) - sata, satz, suna, sunz = optional_datasets # get the dask array underneath sata = sata.data diff --git a/satpy/resample.py b/satpy/resample.py index d9f19073ee..e661ea86ff 100644 --- a/satpy/resample.py +++ b/satpy/resample.py @@ -193,7 +193,7 @@ def hash_dict(the_dict, the_hash=None): """Calculate a hash for a dictionary.""" if the_hash is None: - the_hash = hashlib.sha1() + the_hash = hashlib.sha1() # nosec the_hash.update(json.dumps(the_dict, sort_keys=True).encode('utf-8')) return the_hash @@ -991,7 +991,9 @@ def aggregate(d, y_size, x_size): new_chunks = (tuple(int(x / y_size) for x in d.chunks[0]), tuple(int(x / x_size) for x in d.chunks[1])) - return da.core.map_blocks(_mean, d, y_size, x_size, dtype=d.dtype, chunks=new_chunks) + return da.core.map_blocks(_mean, d, y_size, x_size, + meta=np.array((), dtype=d.dtype), + dtype=d.dtype, chunks=new_chunks) @classmethod def expand_reduce(cls, d_arr, repeats): From c3303e158c2086564bb2fb16514364fecf90cea7 Mon Sep 17 00:00:00 2001 From: David Hoese Date: Fri, 3 Dec 2021 06:46:33 -0600 Subject: [PATCH 6/6] Fix incorrect units in cos sza docstring --- satpy/modifiers/angles.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/satpy/modifiers/angles.py b/satpy/modifiers/angles.py index 3f57e070a2..62af9726ab 100644 --- a/satpy/modifiers/angles.py +++ b/satpy/modifiers/angles.py @@ -268,7 +268,7 @@ def get_satellite_zenith_angle(data_arr: xr.DataArray) -> xr.DataArray: Note that this function can benefit from the ``satpy.config`` parameters :ref:`cache_lonlats ` and :ref:`cache_sensor_angles ` - being set to ``True``. + being set to ``True``. Values are in degrees. """ satz = _get_sensor_angles(data_arr)[1] @@ -279,7 +279,7 @@ def get_cos_sza(data_arr: xr.DataArray) -> xr.DataArray: """Generate the cosine of the solar zenith angle for the provided data. Returns: - DataArray with the same shape as ``data_arr``. Units are radians. + DataArray with the same shape as ``data_arr``. """ lons, lats = _get_valid_lonlats(data_arr.attrs["area"], data_arr.chunks)