Adding a new dtype

[Tom Aldcroft]

By odd coincidence it looks like I got involved in a numpy project to
create a package of extended dtypes. Check the "Quaternion data type"
thread on numpy-discussion. In any case we now have an concrete
example of how to build a new dtype [1]. Not trivial, but not
impossible, and we can more properly evaluate the idea of doing our
datetime object as a dtype.

For now I get the impression there is no user documentation on dtypes
since the quaternion dtype was written at a sprint with a lot of help
from one of the core developers.

- Tom

[1] https://github.com/taldcroft/numpy_quaternion

[Wolfgang Kerzendorf]

Hey Tom,

That is great news. So who was the core developer, who helped you guys? Maybe he'll help us as well ;-)

Cheers
Wolfgang

[Wolfgang Kerzendorf]

Ah by the way I was not only thinking for time but also for celestial coordinates.
But let's start with time.

Cheers,

Wolfgang
On 2012-05-06, at 8:11 AM, Tom Aldcroft wrote:

[Tom Aldcroft]

The core dev who initially helped is Mark Weibe. The quaternion code
and numpy dtype glue is pure C. It isn't entirely clear to me that we
can achieve the sort of flexibility and simplicity that we'd like for
astropy.time in taking the dtype path, but now at least we have an
example that will help making that evaluation.

- Tom

[Wolfgang Kerzendorf]

Hey Tom,

I think we should think a bit more about this dtype thing. Do you want to take this off the list and get back when we have something?

Cheers
Wolfgang

[Erik Bray]

That's too bad that the lack of documentation made things so
difficult. But the implementation itself is, in the end, simple enough
(for a quaternion type anyways). Either way, it makes for a great
example to emulate, should it come to that.

Erik B.

[Erik Tollerud]

I would push back for making coordinates a dtype. It is very
important that coordinates code be fairly readable and usable as a
starting point, given all the affiliated packages out there that
already use various types of coordinate representations. Going down
the dtype road would definitely cut out a lot of the astronomy
developer base that knows python well but is not as C-fluent. There
are also a variety of devil-in-the-detail sorts of problems that
spring to my mind when considering how to do coordinate conversions
using a framework like this. We certainly need to have a
vectorized/array-based coordinate scheme, but there are ways to do
this without necessitating that the implementation be in C.

On the other hand, I do see the rational for this in the case of time
- that's a bit more fixed of an implementation. But I think it's
important we get a time interface done sooner, and doing it the dtype
way will quite likely take much longer...

In the immortal words of Donald Knuth: premature optimization is the
root of all evil!

--
Erik Tollerud
Graduate Student
Center For Cosmology
Department of Physics and Astronomy
2140 Frederick Reines Hall
University of California, Irvine
Cell: (651)307-9409
etol...@uci.edu
http://ps.uci.edu/~etolleru

[Tom Aldcroft]

I agree with Erik on all points below. I've now spent a few hours
messing around with a couple of different dtype implementations in
https://github.com/numpy/numpy-dtypes. There is a huge amount of
overhead and complexity that comes with a dtype. I also discovered
the answer to a question that was nagging me, which is whether custom
attributes in a dtype object get propagated to an array of that dtype.
Example:

>>> import numpy as np
>>> import quaternion
>>> q1 = np.quaternion(1,2,3,4)
>>> q1.w # quaternion components are (w, x, y, z)
1.0
>>> arr = np.array([q1])
>>> arr.dtype
dtype(quaternion)
>>> arr.w
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
AttributeError: 'numpy.ndarray' object has no attribute 'w'

I guess I should not have expected miracles, but I wanted arr.w to
return an array of the "w" attributes. Without that feature then an
array of astropy.Time dtype is less special because you can't
efficiently do vectorized format conversions (e.g. internal => JD)
without resorting to Python wrapper code. So why bother for our case,
where converting representations is most of the battle.

Now *maybe* the quaternion class just isn't doing something right, but
it's development was guided by a core developer (Mark Weibe) at a
sprint, so maybe that's just the way it is. If you look at datetime64
the way to change representation is basically by making a new array
view with a dtype syntax like dtype="datetime_astro[JD]". This won't
fly for our needs, even for datetime. For coords I think it would be
insane.

+MANY on getting things rolling on a straightforward vectorized
solution that we can code today.

- Tom

[Tom Aldcroft]

And just to be clear, by vectorized solution I mean a container class
with an internal array of the float64 or int64/float64 pairs that are
the internal time representation. With just a little effort such a
Time object (or rather it's attributes) can look like either a scalar
or an array with arbitrary shape, depending on how it was initialized.

I meant quaternion dtype.

[Erik Bray]

On Mon, May 7, 2012 at 10:14 PM, Tom Aldcroft
<aldc...@head.cfa.harvard.edu> wrote:

That *should* be possible. If you have an array of complex numbers
then arr.real and arr.imag return properly strided views of the
original array. The quaternion implementation just isn't going
through whatever loops are necessary for that.

But I agree that for now it's probably a moot point.

Erik B.