current community

your communities

more stack exchange communities

Stack Exchange
tour help
careers 2.0
Take the 2-minute tour ×
Stack Overflow is a question and answer site for professional and enthusiast programmers. It's 100% free, no registration required.
up vote 1 down vote favorite
1

I have a large 4D data set and need to create a smaller 4D array from it. I am fairly new to python and am use to IDL or matlab. I read in my values then using the where function I find the index numbers I need for each dimension from smaller 1D arrays. I am trying to create a new array from these index numbers but I keep getting the shape mismatch error (can not be broacast to a single shape.

import numpy as n
import matplotlib.pyplot as plt
import Scientific.IO.NetCDF as S

file=S.NetCDFFile('wspd.mon.mean.nc',mode='r') #Opening File
Lat=file.variables['lat'].getValue()     # Reading in the latitude variables, 73
Lon=file.variables['lon'].getValue()     # Reading in the longitude variables, 144
Level=file.variables['level'].getValue() # Reading in the levels, 17 of them
Defaulttime=file.variables['time'].getValue()   # Reading in the time, hrs since 1-1-1
Defaultwindspeed=file.variables['wspd'].getValue()     # Reading in the windspeed(time, level, lat, lon)

Time=n.arange(len(Defaulttime))/12.+1948  #Creates time array into readable years with 12 months
goodtime=n.where((Time>=1948)&(Time<2013)) #Creates a time array for the years that I want, 1948-2012, since 2013 only has until October, I will not be using that data.
goodlat=n.where((Lat>=35)&(Lat<=50))  #Latitudes where the rockies and plains are in the US
plainslon=n.where((Lon>=275)&(Lon<=285))

Windspeedsplains=Defaultwindspeed[goodtime,:,goodlat,plainslon]

The error below is generated by the line above (last line of code).

>>>ValueError: shape mismatch: objects cannot be broadcast to a single shape
share|improve this question
1  
i doubt anyone will be able to help you with a question that looks like this. The least you could do is figure out which line is actually failing –  Hammer Nov 27 '13 at 1:42
    
Is Defaultwindspeed actually 4D? Use Defaultwindspeed.shape or Defaultwindspeed.ndim to check... –  atomh33ls Nov 27 '13 at 11:40
    
@atomh33ls it must be, otherwise the error would be too many indices. I think the problem is in the shapes of goodtime, goodlat, and plainslon. @Cwilliams, what do you get if you print goodtime.shape, goodlat.shape, plainslon.shape? –  askewchan Nov 27 '13 at 14:03

1 Answer 1

up vote 0 down vote accepted

What is happening is that the lengths of each of your index (where) arrays is going to be different, and so the shape of your output array is ambiguous, hence the error. To force the broadcasting to the proper shape, you'll have to reshape your arrays to broadcast to a proper shape, something like this:

Windspeedsplains = Defaultwindspeed[goodtime[0][:, None, None, None],:,goodlat[0][:,None],plainslon[0]]

The [0] are there becuase np.where(a) returns a tuple of length a.ndim with each element of the tuple being an array of indices that fit your condition. I'm assuming all of your boolean arrays are 1d, so all the where outputs are going to be tuples of length 1, so we just want the one array from it, hence [0].

After we get the array, we want to reshape it to match the shape that you want your output array to have. Presumable, the shape your output should have is (goodtime.size, Defaultwindspeed.shape[1], goodlat.size, plainslon.size), so you must make each of your index arrays have the shape that matches the axis along which the output array should vary for that variable. For example, for goodtime, you want Windspeedplains to vary with respect to time along axis 0 of 4 axes. So, goodtime itself must also vary only along axis 0 of four axes, so you force the index array to have shape (N, 1, 1, 1) which is what [:, None, None, None] does.

So, you could make the above line more readable with:

goodtime = n.where((Time>=1948)&(Time<2013))[0][:, None, None, None]
goodlat = n.where((Lat>=35)&(Lat<=50))[0][:, None]
plainslon = n.where((Lon>=275)&(Lon<=285))[0]

Windspeedsplains=Defaultwindspeed[goodtime, :, goodlat, plainslon]

Or actually, since you can index with the boolean arrays directly:

goodtime = ((Time>=1948)&(Time<2013))[:, None, None, None]
goodlat = ((Lat>=35)&(Lat<=50))[:, None]
plainslon = ((Lon>=275)&(Lon<=285))

Windspeedsplains=Defaultwindspeed[goodtime, :, goodlat, plainslon]

Here is a somewhat simpler example:

In [52]: a = np.arange(3*3*3).reshape(3,3,3)

In [53]: a
Out[53]: 
array([[[ 0,  1,  2],
        [ 3,  4,  5],
        [ 6,  7,  8]],

       [[ 9, 10, 11],
        [12, 13, 14],
        [15, 16, 17]],

       [[18, 19, 20],
        [21, 22, 23],
        [24, 25, 26]]])

In [54]: mask0 = np.where(a[:,0,0] >= 9)

In [55]: mask0
Out[55]: (array([1, 2]),)   # <-- this is the length 1 tuple I was talking about. we want the array inside.

In [56]: mask1 = np.where(a[0,:,0]%2 == 0)

In [57]: mask1
Out[57]: (array([0, 2]),)

In [62]: mask2 = np.where(a[0,0,:] < 1)

In [63]: mask2
Out[63]: (array([0]),)

In [67]: b = a[mask0[0][:, None, None], mask1[0][:, None], mask2[0]]

In [68]: b
Out[68]: 
array([[[ 9],
        [15]],

       [[18],
        [24]]])

In [69]: b.shape
Out[69]: (2, 2, 1)
share|improve this answer

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.