ajdawson / windspharm

A Python library for spherical harmonic computations on vector winds.
http://ajdawson.github.io/windspharm
MIT License
82 stars 36 forks source link

'ValueError: equally-spaced latitudes are invalid (non-global?)' when calculating VectorWind(u, v) #77

Closed Corinne77 closed 7 years ago

Corinne77 commented 7 years ago

I am trying to compute the vorticity according to the 'recipe'

# Read u and v wind components from file.
u = iris.load_cube('uwind_file.nc')
v = iris.load_cube('vwind_file.nc')

# Create an instance of the VectorWind class to do the computations.
w = VectorWind(u, v)
#Call methods to compute streamfunction and relative vorticity.
psi = w.streamfunction()
xi = w.vorticity()

with the difference that I need to concatenate netcdf files (I use CMIP5 data) to cover large time periods.

I get the error message:

w = VectorWind(u,v)
File "/usr/lib/python2.7/site-packages/windspharm/iris.py", line 73, in __init__
 gridtype = self._gridtype(lat.points)
File "/usr/lib/python2.7/site-packages/windspharm/iris.py", line 127, in _gridtype
raise ValueError('equally-spaced latitudes are '
ValueError: equally-spaced latitudes are invalid (non-global?)

I guess that CMIP5 data cover the whole Earth. Can anyone help? Thank you

ajdawson commented 7 years ago

I think this because the IPSL-CM5A-LR model has a grid that includes the poles but has an even number of points in the latitude direction. The underlying software assumes that girds with an even number of points do not include the poles (and conversely grids with an odd number of latitude points must include the poles).

Regridding your data to a suitable grid may be an easy option since you are using iris.

Corinne77 commented 7 years ago

Thank you! I tried regridding using: u.regrid(v, iris.analysis.Linear())

according to the example rotated_air_temp = global_air_temp.regrid(rotated_psl, iris.analysis.Linear()) in http://scitools.org.uk/iris/docs/latest/userguide/interpolation_and_regridding.html,

but I am not sure this is correct. The code seems to take ages or to be stuck doing this...

ajdawson commented 7 years ago

It may take a long time to regrid if you have a lot of data. You may want to look into caching the interpolator.

Corinne77 commented 7 years ago

Thank you. I tried interpolating..

        delta_latitude = 180/96.0
        sample_points = [('longitude', u.coord('longitude').points),('latitude',  np.linspace(90 - 0.5 *  
                                     delta_latitude,-90 + 0.5 * delta_latitude,96))]
        result = u.interpolate(sample_points, iris.analysis.Linear())
        print(result.summary(shorten=True))

which gives me as result eastward_wind / (m s-1) (time: 13140; latitude: 96; longitude: 96) northward_wind / (m s-1) (time: 13140; latitude: 96; longitude: 96)

but now I get the following error message when I want to compute w = VectorWind(u,v) (for u = result): raise ValueError('u and v cannot contain missing values') ValueError: u and v cannot contain missing values

How should I proceed?

ajdawson commented 7 years ago

I've taken a look at the data on JASMIN, it looks like the horizontal interpolation step is successful (and necessary), but the problem is that this data set really does contain missing values. Try plotting the first level at the first time, you'll see what I mean. Windspharm cannot handle data with missing values.

The missing values are there because the data are already interpolated onto pressure levels, and the lower pressure levels intersect the surface in regions of high topography. Some model output chooses to extrapolate below ground (e.g., ERA) but in this case that has not been done.

If I leave out the lowest few pressure levels (the ones that intersect the surface) you can successfully create the VectorWind instance. If you really want the vorticity on the lower levels too, you'll either need to interpolate to fill in the gaps yourself or use a finite difference method to compute the vorticity.