satpy.modifiers.angles module
Utilties for getting various angles for a dataset..
- class satpy.modifiers.angles.ZarrCacheHelper(func: ~typing.Callable, cache_config_key: str, uncacheable_arg_types=(<class 'pyresample.geometry.SwathDefinition'>, <class 'xarray.core.dataarray.DataArray'>, <class 'dask.array.core.Array'>), sanitize_args_func: ~typing.Callable | None = None, cache_version: int = 1)[source]
Bases:
object
Helper for caching function results to on-disk zarr arrays.
It is recommended to use this class through the
cache_to_zarr_if()
decorator rather than using it directly.Currently the cache does not perform any limiting or removal of cache content. That is left up to the user to manage. Caching is based on arguments passed to the decorated function but will only be performed if the arguments are of a certain type (see
uncacheable_arg_types
). The cache value to use is purely based on the hash value of all of the provided arguments along with the “cache version” (see below).Note that the zarr format requires regular chunking of data. That is, chunks must be all the same size per dimension except for the last chunk. To work around this limitation, this class will determine a good regular chunking based on the existing chunking scheme, rechunk the input arguments, and then rechunk the results before returning them to the user. This rechunking is only done if caching is enabled.
- Parameters:
func – Function that will be called to generate the value to cache.
cache_config_key – Name of the boolean
satpy.config
parameter to use to determine if caching should be done.uncacheable_arg_types – Types that if present in the passed arguments should trigger caching to not happen. By default this includes
SwathDefinition
,xr.DataArray
, andda.Array
objects.sanitize_args_func – Optional function to call to sanitize provided arguments before they are considered for caching. This can be used to make arguments more “cacheable” by replacing them with similar values that will result in more cache hits. Note that the sanitized arguments are only passed to the underlying function if caching will be performed, otherwise the original arguments are passed.
cache_version – Version number used to distinguish one version of a decorated function from future versions.
Notes
Caching only supports dask array values.
This helper allows for an additional
cache_dir
parameter to override the use of thesatpy.config
cache_dir
parameter.
Examples
To use through the
cache_to_zarr_if()
decorator:@cache_to_zarr_if("cache_my_stuff") def generate_my_stuff(area_def: AreaDefinition, some_factor: int) -> da.Array: # Generate return my_dask_arr
To use the decorated function:
with satpy.config.set(cache_my_stuff=True): my_stuff = generate_my_stuff(area_def, 5)
Hold on to provided arguments for future use.
- cache_clear(cache_dir: str | None = None)[source]
Remove all on-disk files associated with this function.
Intended to mimic the
functools.cache()
behavior.
- satpy.modifiers.angles._chunks_are_irregular(chunks_tuple: tuple) bool [source]
Determine if an array is irregularly chunked.
Zarr does not support saving data in irregular chunks. Regular chunking is when all chunks are the same size (except for the last one).
- satpy.modifiers.angles._dim_index_with_default(dims: tuple, dim_name: str, default: int) int [source]
- satpy.modifiers.angles._get_output_chunks_from_func_arguments(args)[source]
Determine what the desired output chunks are.
It is assumed a tuple of tuples of integers is defining chunk sizes. If a tuple like this is not found then arguments are checked for array-like objects with a
.chunks
attribute.
- satpy.modifiers.angles._get_sensor_angles(data_arr: DataArray) tuple[DataArray, DataArray] [source]
- satpy.modifiers.angles._get_sensor_angles_ndarray(lons, lats, start_time, sat_lon, sat_lat, sat_alt) ndarray [source]
- satpy.modifiers.angles._get_sun_azimuth_ndarray(lons: ndarray, lats: ndarray, start_time: datetime) ndarray [source]
- satpy.modifiers.angles._hash_args(*args, unhashable_types=(<class 'pyresample.geometry.SwathDefinition'>, <class 'xarray.core.dataarray.DataArray'>, <class 'dask.array.core.Array'>))[source]
- satpy.modifiers.angles._regular_chunks_from_irregular_chunks(old_chunks: tuple[tuple[int, ...], ...]) tuple[tuple[int, ...], ...] [source]
- satpy.modifiers.angles._sunzen_corr_cos_ndarray(data: ndarray, cos_zen: ndarray, limit: float, max_sza: float | None) ndarray [source]
- satpy.modifiers.angles._sunzen_reduction_ndarray(data: ndarray, sunz: ndarray, limit: float, max_sza: float, strength: float) ndarray [source]
- satpy.modifiers.angles.cache_to_zarr_if(cache_config_key: str, uncacheable_arg_types=(<class 'pyresample.geometry.SwathDefinition'>, <class 'xarray.core.dataarray.DataArray'>, <class 'dask.array.core.Array'>), sanitize_args_func: ~typing.Callable | None = None) Callable [source]
Decorate a function and cache the results as a zarr array on disk.
This only happens if the
satpy.config
boolean value for the provided key isTrue
as well as some other conditions. SeeZarrCacheHelper
for more information. Most importantly, this decorator does not limit how many items can be cached and does not clear out old entries. It is up to the user to manage the size of the cache.
- satpy.modifiers.angles.compute_relative_azimuth(sat_azi: DataArray, sun_azi: DataArray) DataArray [source]
Compute the relative azimuth angle.
- Parameters:
sat_azi – DataArray for the satellite azimuth angles, typically in 0-360 degree range.
sun_azi – DataArray for the solar azimuth angles, should be in same range as sat_azi.
- Returns:
A DataArray containing the relative azimuth angle in the 0-180 degree range.
NOTE: Relative azimuth is defined such that: Relative azimuth is 0 when sun and satellite are aligned on one side of a pixel (back scatter). Relative azimuth is 180 when sun and satellite are directly opposite each other (forward scatter).
- satpy.modifiers.angles.get_angles(data_arr: DataArray) tuple[DataArray, DataArray, DataArray, DataArray] [source]
Get sun and satellite azimuth and zenith angles.
Note that this function can benefit from the
satpy.config
parameters cache_lonlats and cache_sensor_angles being set toTrue
.- Parameters:
data_arr – DataArray to get angles for. Information extracted from this object are
.attrs["area"]
,``.attrs[“start_time”]``, and.attrs["orbital_parameters"]
. Seesatpy.utils.get_satpos()
and Metadata for more information. Additionally, the dask array chunk size is used when generating new arrays. The actual data of the object is not used.- Returns:
Four DataArrays representing sensor azimuth angle, sensor zenith angle, solar azimuth angle, and solar zenith angle. All values are in degrees. Sensor angles are provided in the [0, 360] degree range. Solar angles are provided in the [-180, 180] degree range.
- satpy.modifiers.angles.get_cos_sza(data_arr: DataArray) DataArray [source]
Generate the cosine of the solar zenith angle for the provided data.
- Returns:
DataArray with the same shape as
data_arr
.
- satpy.modifiers.angles.get_satellite_zenith_angle(data_arr: DataArray) DataArray [source]
Generate satellite zenith angle for the provided data.
Note that this function can benefit from the
satpy.config
parameters cache_lonlats and cache_sensor_angles being set toTrue
. Values are in degrees.
- satpy.modifiers.angles.sunzen_corr_cos(data: Array, cos_zen: Array, limit: float = 88.0, max_sza: float | None = 95.0) Array [source]
Perform Sun zenith angle correction.
The correction is based on the provided cosine of the zenith angle (
cos_zen
). The correction is limited tolimit
degrees (default: 88.0 degrees). For larger zenith angles, the correction is the same as at thelimit
ifmax_sza
is None. The default behavior is to gradually reduce the correction pastlimit
degrees up tomax_sza
where the correction becomes 0. Bothdata
andcos_zen
should be 2D arrays of the same shape.