If your Kernel class has a predictable amount of member data, then you could define a dtype for it instead of a class. e.g. if it's parameterized by 9 floats and an int, you could do
kerneldt = np.dtype([('myintname', np.int32), ('myfloats', np.float64, 9)])
arr = np.empty(dims, dtype=kerneldt)
You'll have to do some coercion to turn them into objects of class Kernel every time you want to manipulate methods of a single kernel but that's one way to store the actual data in a NumPy array. If you want to only store a reference, then the object dtype is the best you can do without subclassing ndarray.
Answer from dwf on Stack OverflowIf your Kernel class has a predictable amount of member data, then you could define a dtype for it instead of a class. e.g. if it's parameterized by 9 floats and an int, you could do
kerneldt = np.dtype([('myintname', np.int32), ('myfloats', np.float64, 9)])
arr = np.empty(dims, dtype=kerneldt)
You'll have to do some coercion to turn them into objects of class Kernel every time you want to manipulate methods of a single kernel but that's one way to store the actual data in a NumPy array. If you want to only store a reference, then the object dtype is the best you can do without subclassing ndarray.
It has to be a Numpy scalar type:
http://docs.scipy.org/doc/numpy/reference/arrays.scalars.html#arrays-scalars-built-in
or a subclass of ndarray:
http://docs.scipy.org/doc/numpy/reference/generated/numpy.ndarray.html#numpy.ndarray
All 3 requirements conflict with a view.
Ignoring the header field requires selecting the other fields. Selecting a single field is clearly a view, but the state of multiple fields is in flux. When I try anything besides simply viewing the values I get a warning:
In [497]: dt=np.dtype('U10,f,f,f,f')
In [498]: x=np.zeros((5,),dt)
In [505]: x[['f1','f3']].__array_interface__
/usr/bin/ipython3:1: FutureWarning: Numpy has detected that you (may be) writing to an array returned
by numpy.diagonal or by selecting multiple fields in a record
array. This code will likely break in a future numpy release --
see numpy.diagonal or arrays.indexing reference docs for details.
The quick fix is to make an explicit copy (e.g., do
arr.diagonal().copy() or arr[['f0','f1']].copy()).
Remember, the data is layed out element by element, with the dtype tuple values in compact blocks - essentially a compact version of the display. Ignoring the header requires skipping that set of bytes. view can handle skips produced by strides, but not these dtype field skips.
In [533]: x
Out[533]:
array([('header', 0.0, 5.0, 1.0, 10.0), ('header', 1.0, 4.0, 1.0, 10.0),
('header', 2.0, 3.0, 1.0, 10.0), ('header', 3.0, 2.0, 1.0, 10.0),
('header', 4.0, 1.0, 1.0, 10.0)],
dtype=[('f0', '<U10'), ('f1', '<f4'), ('f2', '<f4'), ('f3', '<f4'), ('f4', '<f4')])
To explore reordering the complex fields, lets try a 2d array:
In [509]: y=np.arange(10.).reshape(5,2) # 2 column float
In [510]: y.view(complex) # can be viewed as complex
Out[510]:
array([[ 0.+1.j],
[ 2.+3.j],
[ 4.+5.j],
[ 6.+7.j],
[ 8.+9.j]])
In [511]: y[:,::-1].view(complex)
...
ValueError: new type not compatible with array.
To switch the real/imaginay columns I have to make a copy. complex requires that the 2 floats be contiguous and in order.
In [512]: y[:,::-1].copy().view(complex)
Out[512]:
array([[ 1.+0.j],
[ 3.+2.j],
[ 5.+4.j],
[ 7.+6.j],
[ 9.+8.j]])
float32 to float64 is clearly not a view change. One uses 4 bytes per number, the other 8. You can't 'view' 4 as 8 without copying.
@hpaulj is absolutely correct that this conflicts with a view.
However, you may be asking the wrong question.
numpy can certainly do what you're wanting to do but you'll need to make a temporary copy in memory.
Overall, you're probably better served by rethinking the "read the entire file into memory and then view it" approach. Instead seek past (or read in) the header, then read in the data portion with fromfile. After than, it's relatively straightforward to manipulate things into what you want, as long as you don't mind making a copy to go from float32's to float64's.
To start out with, let's generate a file similar to yours:
import numpy as np
reals = np.arange(100).astype(np.float32)
imag = -9999.0 * np.ones(100).astype(np.float32)
data = np.empty(reals.size + imag.size, dtype=np.float32)
data[::2], data[1::2] = imag, reals
with open('temp.dat', 'wb') as outfile:
# Write a 1Kb header (of literal "x"'s, in this case)
outfile.write(1024 * 'x')
outfile.write(data)
Now we'll read it in.
The key to ignoring the header is to seek past it before reading the data in with fromfile.
Then, we can de-interleave the data and convert to 64-bit floats at the same time.
Finally, you can then view the resulting 2xN-length float64 array as an N-length complex128 array. (Note: complex128 is the 64-bit version of a complex number. complex64 is the 32-bit version.)
For example:
import numpy as np
with open('temp.dat', 'rb') as infile:
# Seek past header
infile.seek(1024)
# Read in rest of file as float32's
data = np.fromfile(infile, dtype=np.float32)
result = np.empty(data.size, np.float64)
# De-interleave imag & real back into expected real & imag, converting to 64-bit
result[::2], result[1::2] = data[1::2], data[::2]
# View the result as complex128's (i.e. 64-bit complex numbers)
result = result.view(np.complex128)
Numpy arrays are most suitable for data types with fixed size. If the objects in the array are not fixed size (such as your MultiEvent) the operations can become much slower.
I would recommend you to store all of the survival times in a 1d linear record array with 3 fields: event_id, time, period. Each event can appear mutliple times in the array:
>>> import numpy as np
>>> rawdata = [(1, 0.4, 4), (1, 0.6, 6), (2,2.6, 6)]
>>> npdata = np.rec.fromrecords(rawdata, names='event_id,time,period')
>>> print npdata
[(1, 0.40000000000000002, 4) (1, 0.59999999999999998, 6) (2, 2.6000000000000001, 6)]
To get data for a specific index you could use fancy indexing:
>>> eventdata = npdata[npdata.event_id==1]
>>> print eventdata
[(1, 0.40000000000000002, 4) (1, 0.59999999999999998, 6)]
The advantage of this approach is that you can easily intergrate it with your ndarray-based functions. You can also access this arrays from cython as described in the manual:
cdef packed struct Event:
np.int32_t event_id
np.float64_t time
np.float64_6 period
def f():
cdef np.ndarray[Event] b = np.zeros(10,
dtype=np.dtype([('event_id', np.int32),
('time', np.float64),
('period', np.float64)]))
<...>
I apologise for not answering the question directly, but I've had similar problems before, and if I understand correctly, the real problem you're now having is that you have variable-length data, which is really, really not one of the strengths of numpy, and is the reason you're running into performance issues. Unless you know in advance the maximum number of entries for a multievent, you'll have problems, and even then you'll be wasting loads of memory/disk space filled with zeros for those events that aren't multi events.
You have data points with more than one field, some of which are related to other fields, and some of which need to be identified in groups. This hints strongly that you should consider a database of some form for storing this information, for performance, memory, space-on-disk and sanity reasons.
It will be much easier for a person new to your code to understand a simple database schema than a complicated, hacked-on-numpy structure that will be frustratingly slow and bloated. SQL queries are quick and easy to write in comparison.
I would suggest based on my understanding of your explanation having Event and MultiEvent tables, where each Event entry has a foreign key into the MultiEvent table where relevant.