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.

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

MV2.repeat unexpected behavior #423

Open
pochedls opened this issue Oct 28, 2020 · 0 comments
Open

MV2.repeat unexpected behavior #423

pochedls opened this issue Oct 28, 2020 · 0 comments

Comments

@pochedls
Copy link

Describe the bug
Behavior of MV2.repeat is unexpected. If I have a transientVariable a = [0, 1, 2, 3], and apply MV2.repeat(a, 2) I get: [0, 0, 1, 1, 2, 2, 3, 3] whereas I would expect [0, 1, 2, 3, 0, 1, 2, 3].

Is MV2.repeat doing the right thing?

To Reproduce

Simple example

[In]: import cdms2
[In]: import MV2
[In]: a = cdms2.createVariable([0, 1, 2, 3])
[In]: MV2.repeat(a, 2)
[Out]: masked_array(data=[0, 0, 1, 1, 2, 2, 3, 3],
             mask=False,
       fill_value=999999)

Applied Example: Remove seasonal cycle from each gridpoint

import cdms2
import MV2
import numpy as np
import matplotlib.pyplot as plt

fn = '/p/user_pub/work/pochedley1/ssmi/ssmi_vapor_1Deg.nc3.nc'
f = cdms2.open(fn)
prw = f('prw')
f.close()
time = prw.getTime()
decimalTime = np.arange(1988, 2020, 1/12.)
lat = prw.getLatitude()
lon = prw.getLongitude()
nyrs = int(len(time)/12)
tropicalInds = np.where(np.abs(lat) < 20)[0]

# MV2 to subtract seasonal cycle
prwClimatology = MV2.reshape(prw, (nyrs, 12, len(lat), len(lon)))  # reshape to (year, month, lat, lon)
prwClimatology = MV2.average(prwClimatology, axis=0)  # average over all years
prwAnom = prw - MV2.repeat(prwClimatology, nyrs, axis=0)  # substract repeated climatology
# simple tropical average (for plotting)
prwAnomTrop = MV2.average(MV2.average(MV2.take(prwAnom, tropicalInds, axis=1), axis=2), axis=1)

# numpy to subtract seasonal cycle
prw_n = np.array(prw)  # to numpy array
prw_n[prw_n < 0] = np.nan  # nan missing data
prw_n[prw_n > 1000] = np.nan  # nan missing data
prwClimatology_n = np.reshape(prw_n, (nyrs, 12, len(lat), len(lon)))  # reshape to (year, month, lat, lon)
prwClimatology_n = np.mean(prwClimatology_n, axis=0)  # average over all years
prwAnom_n = prw_n - np.tile(prwClimatology_n, (nyrs, 1, 1))  # substract repeated climatology
# simple tropical average (for plotting)
prwAnomTrop_n = np.nanmean(np.nanmean(np.take(prwAnom_n, tropicalInds, axis=1), axis=2), axis=1)

# plot two series
plt.plot(decimalTime, np.array(prwAnomTrop), label='MV2')
plt.plot(decimalTime, prwAnomTrop_n, label='numpy')
plt.legend()
plt.xlabel('Time')
plt.ylabel('tropical prw')
plt.show()

test

Note: the numpy version of the time series is correct and the expected behavior.

Expected behavior
In the longer example, I expect MV2.repeat to behave like np.tile (repeating the annual cycle so that it can be subtracted from the prw time series in order to produce the anomaly time series).

Desktop (please complete the following information):

  • OS: Linux
  • Version CDAT8.2 / Python 3
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant