Skip to content
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

More support for pop data #240

Merged
merged 5 commits into from
Mar 4, 2021
Merged

More support for pop data #240

Show file tree
Hide file tree
Changes from 3 commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
496c535
More support for pop data
allibco Mar 1, 2021
271af63
Merge branch 'main' into alli_dev
allibco Mar 1, 2021
b95b149
fix double quited strings
allibco Mar 1, 2021
d5e7c7c
removed image file - in by mistake
allibco Mar 3, 2021
fcc646d
fix plotting for single lat/lon point in pop
allibco Mar 3, 2021
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
119 changes: 31 additions & 88 deletions docs/source/notebooks/POPData.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

120 changes: 34 additions & 86 deletions docs/source/notebooks/TutorialNotebook.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

31 changes: 14 additions & 17 deletions ldcpy/metrics.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -50,24 +50,20 @@ def __init__(
lat_coord_name = ds.cf['latitude'].name
self._lat_coord_name = lat_coord_name

# assum lat/lon (andorider: time, level, lat, lon -
# where time and level may be missing)
if lat_dim_name is None:
if len(ds.dims) == 2:
lat_dim_name = ds.dims[0]
elif len(ds.dims) == 3:
lat_dim_name = ds.dims[1]
elif len(ds.dims) > 3:
lat_dim_name = ds.dims[3]
self._lat_dim_name = lat_dim_name
dd = ds.cf['latitude'].dims
ll = len(dd)
if ll == 1:
if lat_dim_name is None:
lat_dim_name = dd[0]
if lon_dim_name is None:
lon_dim_name = ds.cf['longitude'].dims[0]
elif ll == 2:
if lat_dim_name is None:
lat_dim_name = dd[0]
if lon_dim_name is None:
lon_dim_name = dd[1]

if lon_dim_name is None:
if len(ds.dims) == 2:
lon_dim_name = ds.dims[1]
elif len(ds.dims) == 3:
lon_dim_name = ds.dims[2]
elif len(ds.dims) > 3:
lon_dim_name = ds.dims[4]
self._lat_dim_name = lat_dim_name
self._lon_dim_name = lon_dim_name

# vertical dimension?
Expand Down Expand Up @@ -150,6 +146,7 @@ def _con_var(self, dir, dataset) -> xr.DataArray:
join='override',
)
# con_var = xr.ufuncs.square((o_1 - o_2))
# print(o_1-o_2)
con_var = np.square((o_1 - o_2))
return con_var

Expand Down
29 changes: 12 additions & 17 deletions ldcpy/plot.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -142,23 +142,18 @@ def get_metrics(self, da):
da_data.attrs = da.attrs

# lat/lon dim names are different for ocn and atm
# for now, we assume 2-4 dims time, level, lat, lon
# (where lat and long are always present but the others are optional)
if len(da_data.dims) == 2:
lat_dim = da_data.dims[0]
lon_dim = da_data.dims[1]
elif len(da_data.dims) == 3:
lat_dim = da_data.dims[1]
lon_dim = da_data.dims[2]
elif len(da_data.dims) == 4:
lat_dim = da_data.dims[3]
lon_dim = da_data.dims[4]
else:
print('WARNING lat/lon dims may be off - assuming locations')
lat_dim = da_data.dims[3]
lon_dim = da_data.dims[4]

# print(lat_dim, lon_dim)
dd = da_data.cf['latitude'].dims

ll = len(dd)
if ll == 1:
lat_dim = dd[0]
lon_dim = da_data.cf['longitude'].dims[0]
elif ll == 2:
lat_dim = dd[0]
lon_dim = dd[1]
# if ll = 0, then no lat/lon dims

# print("here", lat_dim, lon_dim)

if self._plot_type in ['spatial']:
metrics_da = lm.DatasetMetrics(da_data, ['time'])
Expand Down
54 changes: 45 additions & 9 deletions ldcpy/util.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -381,14 +381,22 @@ def subset_data(
"""
ds_subset = ds

# print(ds)

if lon_coord_name is None:
lon_coord_name = ds.cf['longitude'].name
if lat_coord_name is None:
lat_coord_name = ds.cf['latitude'].name
# print(lat_coord_name, lon_coord_name)

latdim = ds_subset.cf[lon_coord_name].ndim
# need dim names
dd = ds_subset.cf['latitude'].dims
if latdim == 1:
lat_dim_name = dd[0]
lon_dim_name = ds_subset.cf['longitude'].dims[0]
elif latdim == 2:
lat_dim_name = dd[0]
lon_dim_name = dd[1]
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

same as above


if start is not None and end is not None:
ds_subset = ds_subset.isel({time_dim_name: slice(start, end + 1)})

Expand All @@ -409,14 +417,42 @@ def subset_data(
if vertical_dim_name in ds_subset.dims:
ds_subset = ds_subset.isel({vertical_dim_name: lev})

if lat is not None:
ds_subset = ds_subset.sel(**{lat_coord_name: lat, 'method': 'nearest'})
ds_subset = ds_subset.expand_dims(lat_coord_name)
if latdim == 1:

if lat is not None:
ds_subset = ds_subset.sel(**{lat_coord_name: lat, 'method': 'nearest'})
ds_subset = ds_subset.expand_dims(lat_dim_name)

if lon is not None:
ds_subset = ds_subset.sel(**{lon_coord_name: lon + 180, 'method': 'nearest'})
ds_subset = ds_subset.expand_dims(lon_dim_name)

elif latdim == 2:

if lat is not None:
if lon is not None:

# lat is -90 to 90
# lon should be 0- 360
ad_lon = lon
if ad_lon < 0:
ad_lon = ad_lon + 360

mlat = ds_subset[lat_coord_name].compute()
mlon = ds_subset[lon_coord_name].compute()
# euclidean dist for now....
di = np.sqrt(np.square(ad_lon - mlon) + np.square(lat - mlat))
index = np.where(di == np.min(di))
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
index = np.where(di == np.min(di))
index = np.argmin(di)

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

good to know - thank you!

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So I tried this suggestion, but don't understand the error. Using the simple example from below again:

import numpy as np
import pandas as pd
import xarray as xr
mylat = np.array([[1., 2., 3., 4.], [1.1, 2.1, 3.1, 3.9], [1, 1.9, 2.9, 3.9]])
mylon = np.array([[3., 4., 5., 6.], [3.1, 4.1, 5.1, 5.9], [3, 3.9, 4.9, 5.9]])
mytimes = pd.date_range('2000-01-01', periods=5)
pop_data = np.random.uniform(low=0.0, high=100.0, size=(3,4,5))
ds_subset = xr.DataArray(
    pop_data, 
     coords={
        "time": mytimes,
        "TLAT": (("nlat", "nlon"), mylat, {'standard_name': 'latitude', 'units': 'degrees_north'}),
        "TLON": (("nlat", "nlon"), mylon, {'standard_name': 'longitude', 'units': 'degrees_east'})
    },
    dims=['nlat', 'nlon', 'time'] , attrs={'long_name': 'Surface Potential'}
)

Then I do this calculation:

mlat = ds_subset["TLAT"].compute()
mlon = ds_subset["TLON"].compute()
ad_lat = 3.1
ad_lon = 4
di = np.sqrt(np.square(ad_lon - mlon) + np.square(ad_lat - mlat))
index = np.argmin(di)

Then I get this error:
"ValueError: dimensions ('nlat', 'nlon') must have the same length as the number of data dimensions, ndim=0 "
(doing 'index = np.where(di == np.min(di))' does not give an error)

xmin = index[0][0]
ymin = index[1][0]
ds_subset = ds_subset.isel({lat_dim_name: xmin, lon_dim_name: ymin})

ds_subset = ds_subset.expand_dims(lat_dim_name)
ds_subset = ds_subset.expand_dims(lon_dim_name)

if lon is not None:
ds_subset = ds_subset.sel(**{lon_coord_name: lon + 180, 'method': 'nearest'})
ds_subset = ds_subset.expand_dims(lon_coord_name)
# print(ds_subset)

# print(ds_subset)
# I need to have TLAT and TLON coords have dims (nlat, nlon) and I
# can't figure this out....

return ds_subset
128 changes: 65 additions & 63 deletions tests/test_metrics.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -112,9 +112,9 @@ def test_mean_error_02(self):
self.assertTrue(em.mean.all() == 0.0)

def test_dim_names(self):
self.assertTrue(test_spatial_metrics._lat_dim_name == "lat")
self.assertTrue(test_spatial_metrics._lon_dim_name == "lon")
self.assertTrue(test_spatial_metrics._time_dim_name == "time")
self.assertTrue(test_spatial_metrics._lat_dim_name == 'lat')
self.assertTrue(test_spatial_metrics._lon_dim_name == 'lon')
self.assertTrue(test_spatial_metrics._time_dim_name == 'time')

def test_TS_02(self):
import xarray as xr
Expand Down Expand Up @@ -147,11 +147,11 @@ def test_min_val(self):
def test_max_val(self):
self.assertTrue(test_overall_metrics.max_val == 99)

# def test_ns_con_var(self):
# self.assertTrue(test_overall_metrics.ns_con_var == 2500) # is this right?
def test_ns_con_var(self):
self.assertTrue(test_overall_metrics.ns_con_var == 2500) # is this right?

# def test_ew_con_var(self):
# self.assertTrue(test_overall_metrics.ew_con_var == 400) # is this right?
def test_ew_con_var(self):
self.assertTrue(test_overall_metrics.ew_con_var == 400) # is this right?

def test_odds_positive(self):
self.assertTrue(np.isclose(test_overall_metrics.odds_positive, 0.98019802, rtol=1e-09))
Expand Down Expand Up @@ -291,19 +291,19 @@ def test_max_val_spatial(self):
).all()
)

# def test_ns_con_var_spatial(self):
# self.assertTrue(
# (
# test_spatial_metrics.get_metric('ns_con_var')
# == np.array(
# [
# [2500.0, 2500.0, 2500.0, 2500.0, 2500.0],
# [2500.0, 2500.0, 2500.0, 2500.0, 2500.0],
# [2500.0, 2500.0, 2500.0, 2500.0, 2500.0],
# ]
# )
# ).all()
# )
def test_ns_con_var_spatial(self):
self.assertTrue(
(
test_spatial_metrics.get_metric('ns_con_var')
== np.array(
[
[2500.0, 2500.0, 2500.0, 2500.0, 2500.0],
[2500.0, 2500.0, 2500.0, 2500.0, 2500.0],
[2500.0, 2500.0, 2500.0, 2500.0, 2500.0],
]
)
).all()
)

def test_odds_positive_spatial(self):
self.assertTrue(
Expand Down Expand Up @@ -446,20 +446,20 @@ def test_zscore_spatial(self):
).all()
)

# def test_ew_con_var_spatial(self):
# self.assertTrue(
# (
# test_spatial_metrics.get_metric('ew_con_var')
# == np.array(
# [
# [100.0, 100.0, 100.0, 100.0, 1600.0],
# [100.0, 100.0, 100.0, 100.0, 1600.0],
# [100.0, 100.0, 100.0, 100.0, 1600.0],
# [100.0, 100.0, 100.0, 100.0, 1600.0],
# ]
# )
# ).all()
# )
def test_ew_con_var_spatial(self):
self.assertTrue(
(
test_spatial_metrics.get_metric('ew_con_var')
== np.array(
[
[100.0, 100.0, 100.0, 100.0, 1600.0],
[100.0, 100.0, 100.0, 100.0, 1600.0],
[100.0, 100.0, 100.0, 100.0, 1600.0],
[100.0, 100.0, 100.0, 100.0, 1600.0],
]
)
).all()
)

def test_mean_time_series(self):
self.assertTrue(
Expand Down Expand Up @@ -524,27 +524,27 @@ def test_min_val_time_series(self):
).all()
)

# def test_ns_con_var_time_series(self):
# self.assertTrue(
# np.isclose(
# test_time_series_metrics.get_metric('ns_con_var'),
# np.array(
# [
# 2500.0,
# 2500.0,
# 2500.0,
# 2500.0,
# 2500.0,
# 2500.0,
# 2500.0,
# 2500.0,
# 2500.0,
# 2500.0,
# ]
# ),
# rtol=1e-09,
# ).all()
# )
def test_ns_con_var_time_series(self):
self.assertTrue(
np.isclose(
test_time_series_metrics.get_metric('ns_con_var'),
np.array(
[
2500.0,
2500.0,
2500.0,
2500.0,
2500.0,
2500.0,
2500.0,
2500.0,
2500.0,
2500.0,
]
),
rtol=1e-09,
).all()
)

def test_odds_positive_time_series(self):
self.assertTrue(
Expand Down Expand Up @@ -679,14 +679,16 @@ def test_zscore_time_series(self):
).all()
)

# def test_ew_con_var_time_series(self):
# self.assertTrue(
# np.isclose(
# test_time_series_metrics.get_metric('ew_con_var'),
# np.array([400.0, 400.0, 400.0, 400.0, 400.0, 400.0, 400.0, 400.0, 400.0, 400.0]),
# rtol=1e-09,
# ).all()
# )
def test_ew_con_var_time_series(self):
self.assertTrue(
np.isclose(
test_time_series_metrics.get_metric('ew_con_var'),
np.array(
[400.0, 400.0, 400.0, 400.0, 400.0, 400.0, 400.0, 400.0, 400.0, 400.0]
),
rtol=1e-09,
).all()
)

def test_diff_pcc(self):
self.assertTrue(
Expand Down
Loading