Turn a cube into a list of data points and lists of x, y and z

Question

This is a bit of a long shot, but I thought I'd ask here before writing this myself.

I have a 3D cube of data with lon, lat and height coordinates. I want 4 1D vectors of the data at all points, the lon, lat, height and data. This is so that I can then write it to an ASCII file as a list of points with their locations. Doing this for the data is easy with a reshape, but the part that is trickier is turning the coordinates into the right vectors.

Has anyone done this already and have some hints?

score 1 · Accepted Answer · answered Oct 16 '17 at 13:04

1

I'm not 100% sure I've understood the aim correctly, but would something like this do what you want?

lats = []
lons = []
heights = []
data = []
for point_cube in cube.slices_over(['latitude', 'longitude', 'height']):
    lats.append(point_cube.coord('latitude').points[0])
    lons.append(point_cube.coord('longitude').points[0])
    heights.append(point_cube.coord('height').points[0])
    data.append(point_cube.data)

Or for something (almost certainly) more efficient, you could explore using the numpy.meshgrid function to turn your 1-d coord.points arrays into 3-d arrays, which you could then handle in the same way as the data array.

answered Oct 16 '17 at 13:04

RuthC

1,269
1
10
21

Thanks, that would work. I think that meshgrid would be useful as well, good call. – Nick Savage Oct 16 '17 at 13:10
using meshgrid, the problem I am having is in making sure that I will really have the right coordinates for each data point. So the data will simply be given by cube.data.reshape(). However, I don't know if cube data can be relied to be in a specific order (c or fortran) or how to make sure that order of arguments to mesgrid coresponds to this. It is all getting rather complicated... – Nick Savage Oct 23 '17 at 10:06

score 1 · Answer 2 · answered Nov 30 '17 at 14:29

Another option is to use itertools to build the product of all coordinate points and flatten the cube's data array:

points_prod = itertools.product(cube.coord('height').points,
                                cube.coord('latitude').points,
                                cube.coord('longitude').points)
flat_data = cube.data.reshape(-1)

The itertools product is a generator, which we can consume by converting it into a list. We can then index the list of products and the flattened data array. The index i of the list of products will be the coordinate points for the data value i in the flattened data array:

points_prod_list = list(points_prod)
print '{} -- {}'.format(points_prod_list[i], flat_data[i])

Thinking about specific data orderings and ensuring you have the right coordinates for each data point...

So long as you get the order of the dimensions the same as the order that they are printed, this will always link the right coordinate point values to the data value that those point values describe.

This is quite hard to describe but can be demonstrated reasonably simply using the code above. It relates to the fact that the first printed dimension coordinate describes the outermost axis of the cube's data array, and the orders in which itertools.product creates the product of all input combinations, and the order in which NumPy flattens an array. Fundamentally though, Iris itself is relying on these orderings to make sure its coordinate points values describe the right data values!

Turn a cube into a list of data points and lists of x, y and z

2 Answers2