From 54955f73ec994c0f2e4b2387f1d38a2be777a214 Mon Sep 17 00:00:00 2001
From: tomrink <rink@ssec.wisc.edu>
Date: Fri, 25 Oct 2024 11:58:36 -0500
Subject: [PATCH] snapshot...
---
modules/util/resample.py | 44 +++++++++++++++++++++++-----------------
1 file changed, 25 insertions(+), 19 deletions(-)
diff --git a/modules/util/resample.py b/modules/util/resample.py
index 86146634..c25f6e20 100644
--- a/modules/util/resample.py
+++ b/modules/util/resample.py
@@ -1,10 +1,8 @@
import h5py
import numpy as np
from scipy.interpolate import RegularGridInterpolator, griddata
-from cartopy.crs import LambertAzimuthalEqualArea
import cartopy.crs as ccrs
-from netCDF4 import Dataset
-import xarray as xr
+
from util.util import get_grid_values_all
# resample methods:
@@ -72,7 +70,7 @@ def reproject(fld_2d, lat_2d, lon_2d, proj, target_grid=None, grid_spacing=15000
automatically. The target grid is always returned.
:param grid_spacing: distance between the target grid points (in meters)
:param method: resampling method: 'linear', 'nearest', 'cubic'
- :return: reprojected 2D scalar field, the target grid (will be 2D if rotate=True)
+ :return: reprojected 2D scalar field, the target grid
"""
data_xy = proj.transform_points(ccrs.PlateCarree(), lon_2d, lat_2d)[..., :2]
@@ -112,29 +110,37 @@ def bisect_great_circle(lon_a, lat_a, lon_b, lat_b):
return lon_c, lat_c
-def acspo_sst(filename):
+def acspo_sst(filename, stride=2, grid_spacing=2000):
h5f = h5py.File(filename, 'r')
- sst = get_grid_values_all(h5f, 'sea_surface_temperature')
- sst = sst[0, :, :]
- print(sst.shape)
- lons = get_grid_values_all(h5f, 'lon')
- lats = get_grid_values_all(h5f, 'lat')
- print(lons.shape, lats.shape)
- sst = sst[::2, ::2]
- lons = lons[::2, ::2]
- lats = lats[::2, ::2]
+ # operational
+ # sst = get_grid_values_all(h5f, 'sea_surface_temperature')
+ # sst = sst[0, :, :]
+ # lons = get_grid_values_all(h5f, 'lon')
+ # lats = get_grid_values_all(h5f, 'lat')
+
+ # classic
+ sst = get_grid_values_all(h5f, 'sst_regression')
+ sst = sst[:, :]
+ lons = get_grid_values_all(h5f, 'longitude')
+ lats = get_grid_values_all(h5f, 'latitude')
+
+ print('data dims: ', lons.shape, lats.shape)
+
+ sst = sst[::stride, ::stride]
+ lons = lons[::stride, ::stride]
+ lats = lats[::stride, ::stride]
ylen, xlen = lons.shape
- print(ylen, xlen)
+ print('final dims: ', ylen, xlen)
- cen_lon = lons[ylen // 2, xlen // 2]
- cen_lat = lats[ylen // 2, xlen // 2]
- print(cen_lat, cen_lon)
+ cen_lon = lons[ylen // stride, xlen // stride]
+ cen_lat = lats[ylen // stride, xlen // stride]
+ print('center latitude/longitude: ', cen_lat, cen_lon)
proj = get_projection('LambertAzimuthalEqualArea', cen_lat, cen_lon)
- fld_reproj, (y_map_2d, x_map_2d) = reproject(sst, lats, lons, proj, grid_spacing=2000.0)
+ fld_reproj, (y_map_2d, x_map_2d) = reproject(sst, lats, lons, proj, grid_spacing=grid_spacing)
return fld_reproj, proj, y_map_2d, x_map_2d
--
GitLab