Original author: rokroskar (November 29, 2010 15:55:50)

Attempting to run pynbody using python 2.5.4 - seems that the way units are being implemented isn't backwards compatible - are we ok requiring python > 2.6?:

[roskar@miniature pynbody]$ python
Python 2.5.4 (r254:67916, Aug 11 2009, 13:44:50)
[GCC 4.1.2 20071124 (Red Hat 4.1.2-42)] on linux2
Type "help", "copyright", "credits" or "license" for more information.
Traceback (most recent call last):
File "/home/itp/roskar/.pythonrc.py", line 7, in <module>
import pynbody, pynbody.analysis
File "/home/itp/roskar/nbody/python/pynbody/pynbody/init.py", line 1, in <module>
import util, filt, array, decorate, family, snapshot, tipsy, gadget, analysis
File "/home/itp/roskar/nbody/python/pynbody/pynbody/array.py", line 3, in <module>
class SimArray(np.ndarray) :
File "/home/itp/roskar/nbody/python/pynbody/pynbody/array.py", line 13, in SimArray
@units.setter
AttributeError: 'property' object has no attribute 'setter'

Original issue: http://code.google.com/p/pynbody/issues/detail?id=4

Original author: [email protected] (March 08, 2011 03:05:24)

Not sure where I should put this, but I'm just curious how to do something. Perhaps we should have a "general questions" section? Or maybe that is the point of the pynbody google group? Anyway...

If I have an array of arrays, such as s['vel'], and I would like to make a new array representing the vector length of each sub array, how would I do this? I've tried 2 things that I thought would work:

mvl is mean velocity, obtained from 'smooth' program from uwashington. Would like # to calculate mean speed.

s['msp'] = numpy.sqrt(numpy.vdot(s['mvl'],s['mvl']))

AttributeError Traceback (most recent call last)

/2/home/woodsrm/Data/Gas/g15784.new/<ipython console> in <module>()

/home/woodsrm/SimTools/pynbody/pynbody/plot/profile.py in rotation_curve(s, center, r_units, v_units, disk_height, nbins, bin_spacing)
20 pro = profile.Profile(s, nbins = nbins, type=bin_spacing)
21
---> 22 r = pro.r.in_units(r_units)
23 v = pro['v_circ'].in_units(v_units)
24

AttributeError: Profile instance has no attribute 'r'

Also tried:

s['msp'] = [numpy.sqrt(numpy.vdot(i,i)) for i in s['mvl']]

AttributeError Traceback (most recent call last)

/2/home/woodsrm/Data/Gas/g15784.new/<ipython console> in <module>()

/home/woodsrm/SimTools/pynbody/pynbody/snapshot.py in setitem(self, name, item)
128 except TypeError :
129 ndim = 1
--> 130 self._create_array(name, ndim, dtype=item.dtype)
131
132 # Copy in contents if the contents isn't actually pointing to

AttributeError: 'list' object has no attribute 'dtype'

How should I go about doing this? Thanks!

Original issue: http://code.google.com/p/pynbody/issues/detail?id=31

Original author: [email protected] (March 06, 2011 20:41:20)

I don't know when this happened, but the current version of AHFstep renumbers tipsy particles in their particles file to dark, star, gas order instead of the usual gas, dark, star. It seems like we should provide support for either, but I'm not sure how to make the switch.

I will email Steffen and Alex and see what they have to say about this (and advertise pynbody to them as well).

Original issue: http://code.google.com/p/pynbody/issues/detail?id=29

Original author: rokroskar (February 25, 2011 01:46:15)

i tried checking out jonathan's partial loader branch and hg puked all over me - what's the correct order of operations in this case?

Original issue: http://code.google.com/p/pynbody/issues/detail?id=27

Original author: rokroskar (November 18, 2010 14:28:42)

from a UI point of view, it might be nice to be able to specify 'all' as a family. I.e. when making a profile you want to be able to pass a family argument, but if you pass 'all' each routine then needs an if fam is 'all'... etc. I'm not sure how this should be implemented... ?

Original issue: http://code.google.com/p/pynbody/issues/detail?id=3

Original author: rokroskar (March 18, 2011 14:10:17)

should routines in analysis/halo.py that do the centering also return the center of mass velocity? I've added this to the center_of_mass function, but how should velocity be handled in the other ones (i.e. when centering on potential or shrink-sphere or both). Ideally, the center() function would recenter both, the positions and the velocities.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=41

Original author: rokroskar (March 18, 2011 14:25:57)

apparently, the units parser can sometimes switch the order of units, causing the test_unit_parser test to fail:

..........F...

FAIL: units_test.test_units_parser

Traceback (most recent call last):
File "/usr/local/Python/2.6.2/lib/python2.6/site-packages/nose/case.py", line 183, in runTest
self.test(self.arg)
File "/home/itp/roskar/nbody/python/pynbody/nose/units_test.py", line 29, in test_units_parser
assert str(testunit)=="kpc__2/3 a__2/3 s_-2/3 Myr__-2/3"
AssertionError:
-------------------- >> begin captured stdout << ---------------------
kpc__2/3 a__2/3 Myr__-2/3 s_*-2/3

--------------------- >> end captured stdout << ----------------------

Ran 14 tests in 6.284s

FAILED (failures=1)

Original issue: http://code.google.com/p/pynbody/issues/detail?id=42

Original author: [email protected] (January 19, 2011 18:41:45)

e.g. f.gas["smooth"].recalculate() fails because f.gas["smooth"] doesn't know it's a family array

I'll fix this, I'm just leaving this note for myself.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=22

Original author: [email protected] (March 03, 2011 16:37:46)

What steps will reproduce the problem?
In [1]: import pynbody
In [2]: s = pynbody.load('g15784.00256')
In [3]: f = s.halos()

What is the expected output? What do you see instead?

Expected: no output, list of halos

received:

IndexError Traceback (most recent call last)

/2/home/woodsrm/Data/Gas/g15784.new/<ipython console> in <module>()

/home/woodsrm/SimTools/pynbody/pynbody/snapshot.py in halos(self, *args)
143
144 for c in halo._halo_classes :
--> 145 if c._can_load(self) : return c(self, *args)
146
147 raise RuntimeError("No halo catalogue found")

/home/woodsrm/SimTools/pynbody/pynbody/halo.py in init(self, f)
40 HaloCatalogue.init(self)
41 self._load_ahf_halos(glob.glob(f._filename+'_z_halos')[0])
---> 42 self._load_ahf_substructure(glob.glob(f._filename+'_z_substructure')[0])
43 self._particles_loaded = False
44

IndexError: list index out of range

Please use labels and text to provide additional information.

Not sure if this is just something of mine, or if something got broken on a recent update as this used to work.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=28

Original author: rokroskar (December 13, 2010 17:36:34)

I can't figure this one out...

>>> s = pynbody.load('10/12M_hr.01000')
Loading units from /lustre/zbox3b/roskar/isolated_runs/12M_hr/12M_hr.param
>>> s.star['smooth']
SimArray([ NaN, NaN, NaN, ..., NaN, NaN, NaN], 'kpc')
>>> s.star.calculate_smoothing()
>>> s.star['smooth']
SimArray([ 0.00699206, 0.01383823, 0.03779312, ..., 0.04534013,
0.10665303, 0.13607457], 'kpc')

I don't know why the lazy-load isn't working properly, please advise!

Original issue: http://code.google.com/p/pynbody/issues/detail?id=15

Original author: [email protected] (April 02, 2011 14:04:00)

Thanks for commiting g15784 Greg, but it's missing a halo catalogue.

Once you've provided a halo catalogue I can update snapshot_test.py so we don't need MW.512 any more.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=50

Original author: [email protected] (March 23, 2011 18:51:54)

e.g. numpy.where(f['x']<0.5) returns an array of indices with the same units as f['x'].

This is odd, and either a pure numpy array should be returned or a pynbody array with no units.

I will fix this issue, but put this placeholder issue to remind me.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=45

Original author: [email protected] (March 08, 2011 13:46:48)

I got in above my depth when I created the Halo object. It is inherited from IndexedSubSnap, which is inherited down a long chain from SimSnap, so we expect it to act like a SimSnap. And it does in some cases

In [8]: h[1].computable_keys()
Out[8]:
['hetot',
'vrxy',
'ke',
'ofe',
'smooth',
'j',
'vcxy',
'v2',
'feh',
'rxy',
'vr',
'r',
'rho',
'hydrogen']

But then fails in places where it would be really nice for it to work:

===Failure 1===

In [9]: h[1].star['feh']

KeyError Traceback (most recent call last)

/shadow/data2/gsstinson/g15784/bigstarlog/00185/<ipython console> in <module>()

/shadow/data1/gss/home/src/pynbody/pynbody/snapshot.pyc in getitem(self, i)
92 # All available methods of getting this array have failed

 93 

---> 94 raise KeyError, "Unknown array "+i
95
96 elif isinstance(i,slice) :

KeyError: 'Unknown array feh'

===Failure 2===

In [11]: pp.schmidtlaw(h[1])
ERROR: An unexpected error occurred while tokenizing input
The following traceback may be corrupted or invalid
The error message is: ('EOF in multi-line statement', (11, 0))

ERROR: An unexpected error occurred while tokenizing input
The following traceback may be corrupted or invalid
The error message is: ('EOF in multi-line statement', (28, 0))

TypeError Traceback (most recent call last)

/shadow/data2/gsstinson/g15784/bigstarlog/00185/<ipython console> in <module>()

/shadow/data1/gss/home/src/pynbody/pynbody/plot/stars.pyc in schmidtlaw(sim, center, filename, pretime, diskheight, rmax, radial, **kwargs)
37
38 if center :
---> 39 angmom.faceon(sim)
40
41 # select stuff

/shadow/data1/gss/home/src/pynbody/pynbody/analysis/angmom.pyc in faceon(h, *_kwargs)
98 the halo of interest), this routine centers the simulation and
99 rotates it so that the disk lies in the x-y plane. This gives
100 a face-on view for SPH images, for instance."""
101
--> 102 sideon(h, calc_faceon_matrix, *_kwargs)

/shadow/data1/gss/home/src/pynbody/pynbody/analysis/angmom.pyc in sideon(h, vec_to_xform, cen_size, disk_size, cen, vcen, verbose)
52 if verbose :
53 print "Finding halo center..."
---> 54 cen = halo.hybrid_center(h) # or h['pos'][h['phi'].argmin()]
55 if verbose :
56 print "cen=",cen

/shadow/data1/gss/home/src/pynbody/pynbody/analysis/halo.pyc in hybrid_center(sim, r, *_kwargs)
60
61 cen_a = potential_minimum(sim)
---> 62 return shrink_sphere_center(sim[filt.Sphere(r, cen_a)], *_kwargs)
63
64 def center(sim, mode='pot') :

/shadow/data1/gss/home/src/pynbody/pynbody/analysis/halo.pyc in shrink_sphere_center(sim, r, shrink_factor, min_particles, verbose)
17 # use rough estimate for a maximum radius

 18         # results will be insensitive to the exact value chosen

---> 19 r = (sim["x"].max()-sim["x"].min())/2
20
21 while len(x)>min_particles :

/shadow/data1/gss/home/src/pynbody/pynbody/array.pyc in q(s, _y, *_kw)
610 return w(SimArray(s), _y, *_kw)
611 else :
--> 612 return w(SimArray(s), *y)
613
614 q.name = w.name

/shadow/data1/gss/home/src/pynbody/pynbody/array.pyc in new(subtype, data, units, sim, *_kwargs)
143
144 def new(subtype, data, units=None, sim=None, *_kwargs) :
--> 145 new = np.array(data, **kwargs).view(subtype)
146
147 if hasattr(data, 'units') and hasattr(data, 'sim') and units is None and sim is None :

/shadow/data1/gss/home/src/pynbody/pynbody/array.pyc in array(self, _args, *_kwargs)
527
528 def array(self, _args, *_kwargs) :
--> 529 return self.base[self._ptr].array(_args, *_kwargs)
530
531 def _reexpress_index(self, index) :

/shadow/data1/gss/home/src/pynbody/pynbody/array.pyc in getitem(self, item)
536
537 def getitem(self, item) :
--> 538 return self.base[self._reexpress_index(item)]
539
540 def setitem(self, item, to) :

/shadow/data1/gss/home/src/pynbody/pynbody/array.pyc in _reexpress_index(self, index)
533 return [self._ptr[index[0]]]+list(index[1:])
534 else :
--> 535 return self._ptr[index]
536
537 def getitem(self, item) :

TypeError: only integer arrays with one element can be converted to an index

Original issue: http://code.google.com/p/pynbody/issues/detail?id=32

Original author: [email protected] (March 14, 2011 20:16:49)

I wanted to report this just as a heads up as I'm not sure if this behavior is intended.

When using the "faceon" routine, if it is called twice, the 2nd image produced is now rotated counter-clockwise by 90 degrees. Another call will rotate it another 90 degrees, and so on.

Presumably people won't be calling this multiple times, so it's probably not a big deal. However, if someone accidentally calls it twice, it's annoying that the image is now rotated from all of your others.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=37

Original author: rokroskar (December 13, 2010 10:56:52)

this doesn't work:

>>> np.sum(s['mass'])

this does:

>>> s['mass'].sum()

this gives the wrong units:

>>> np.sqrt((s['pos'][:,0:3]2).sum(axis=1))
SimArray([ 150.93175958, 67.12237074, 7.0066101 , ..., 6.47421208,
10.23414716, 10.68530362], 'kpc2')

this gives the correct units:

>>> (s['pos']2).sum(axis=1)(1,2)
SimArray([ 150.93175958, 67.12237074, 7.0066101 , ..., 6.47421208,
10.23414716, 10.68530362], 'kpc')

It would be nice if these sets of operations gave the same results, otherwise this is bound to lead to mayhem...

Original issue: http://code.google.com/p/pynbody/issues/detail?id=11

Original author: [email protected] (March 16, 2011 21:08:45)

This is an odd one... hopefully I haven't just done something dumb.

What steps will reproduce the problem?

In [1]: import pynbody
In [2]: s = pynbody.load('g15784.00168')
In [8]: s.physical_units()

are there any particles in -0.5<x,y<0.5 kpc ?

In [9]: i = numpy.where(s.gas['x'] < 0.5)
In [10]: i = numpy.where(s.gas['x'][i] > -0.5)
In [11]: i = numpy.where(s.gas['y'][i] > -0.5)
In [12]: i = numpy.where(s.gas['y'][i] < 0.5)
In [13]: i
Out[13]:
(SimArray([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 17, 18,
19, 20, 21, 22, 25, 26, 27, 28, 29, 31, 33, 34, 35, 36, 37, 38, 41,
42, 44, 48, 49, 50, 51, 52, 53, 56, 57, 58, 59, 62, 63, 64, 67, 69,
71, 72, 74, 75, 76, 77, 78, 80, 81, 82], 'kpc'),)

yes there are. Now calculate mach number

In [16]: s['mach'] = s['dsp']/s['c']

tried recalculating i again, same as above, reports same thing, so error hasn't occurred yet

In [18]: import pynbody.plot.sph as sph
In [24]: h = s.halos()
In [25]: import pynbody.analysis.angmom as angmom
In [26]: angmom.faceon(h[1])
In [19]: sph.image(s.gas,qty='mach')

weird image, see attached

are there particles in same xy region?

In [30]: i = numpy.where(s.gas['x'] < 0.5)
In [31]: i = numpy.where(s.gas['x'][i] > -0.5)
In [32]: i = numpy.where(s.gas['y'][i] > -0.5)
In [33]: i = numpy.where(s.gas['y'][i] < 0.5)
In [34]: i
Out[34]: (SimArray([], dtype=int64, 'kpc'),)

no?? where did they go??

What is the expected output? What do you see instead?

expected: pretty plot, no change to particles information
instead: plot with missing squares, particles no longer in some regions?

Please use labels and text to provide additional information.

Checking length of gas array before and after shows that nothing has left the array.... but there are no longer particles in this region?? Also checked 'mach' values in this region:
In [35]: s['mach'][i]
Out[35]:
SimArray([ 28.30397303, 26.20634357, 28.75025244, 21.82961108,
21.1740042 , 14.71155915, 25.0539672 , 17.53637813,
26.26491879, 20.9054888 , 15.69044442, 28.83411271,
16.77683446, 14.35983621, 10.18237936, 11.33148468,
18.82266313, 18.10151572, 16.00618036, 16.62856109,
21.4643663 , 24.22581265, 33.78212089, 20.85854108,
22.52107978, 23.43231128, 33.47977775, 14.8898306 ,
13.76156623, 19.17479612, 24.36452344, 28.53118911,
16.52399683, 14.53370011, 18.50623492, 14.93391785,
25.85545107, 14.19155804, 16.05767424, 12.00926287,
13.28488 , 22.13641185, 19.12379424, 13.95318768,
15.66015145, 14.28708515, 9.73865815, 11.72071385,
22.36727074, 11.13602548, 11.08713528, 16.58364836,
11.3453358 , 20.10502023, 11.59132031, 10.05498969,
7.99850095, 12.08812109, 12.03373314, 15.31968273, 15.18570174])

Nothing out of the ordinary! I can provide test files if necessary...

Original issue: http://code.google.com/p/pynbody/issues/detail?id=40

Original author: [email protected] (January 05, 2011 17:02:48)

I just installed python 2.6 from python.org and then scipy0.8 on top of that and got the following, which seems strange:

kemosabe:12M stinson$ python
Python 2.6.1 (r261:67515, Jun 24 2010, 21:47:49)
[GCC 4.2.1 (Apple Inc. build 5646)] on darwin
Type "help", "copyright", "credits" or "license" for more information.
>>> import pynbody
>>> s = pynbody.load('12M.00350')
Loading units from /Users/stinson/12M/12M.param
>>> import pynbody.plot.sph
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/Users/stinson/code/pynbody/pynbody/plot/sph.py", line 3, in <module>
from .. import sph
File "/Users/stinson/code/pynbody/pynbody/sph.py", line 210, in <module>
@snapshot.SimSnap.derived_quantity
NameError: name 'snapshot' is not defined

Original issue: http://code.google.com/p/pynbody/issues/detail?id=17

Original author: [email protected] (January 20, 2011 00:07:04)

>>> import pynbody
>>> s = pynbody.load('g25271.00880')
>>> import pynbody.plot.sph as sphplot
>>> sphplot.image(s.gas,qty='temp')
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/shadow/data1/gss/home/src/pynbody/pynbody/plot/sph.py", line 77, in image
im[np.where(im==0)] = abs(im[np.where(im!=0)]).min()
ValueError: zero-size array to ufunc.reduce without identity

This doesn't seem to happen to the extent of an error with a z=0 output, but produces a blank (all blue) image. I'm still not sure how to set break points.

I'm pretty sure I've got all the code correct now.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=24

Original author: [email protected] (April 01, 2011 08:52:26)

Bug report from Fabio.

On elektra:

/net/astro-quinn/fabio/REPOSITORY/e12Gals/h258.cosmo50cmb.1536gst1MbwKBH/h258.cosmo50cmb.1536gst1MbwKBH.00512

gas density plots attached have a big central 'blob' which is presumably related to the presence of a BH there. However the BH is officially a 'star' particle so how is it ending up causing this disruption in the gas density plots?

Original issue: http://code.google.com/p/pynbody/issues/detail?id=49

Original author: [email protected] (March 14, 2011 16:33:46)

What steps will reproduce the problem?
import pynbody, pynbody.analysis.angmom as angmom
s = pynbody.load()
h = s.halos()
angmom.faceon(h[1])

What is the expected output? What do you see instead?
expected: no output, properly oriented halo

received:

AttributeError Traceback (most recent call last)

/2/home/woodsrm/Data/Gas/g15784.new/<ipython console> in <module>()

/home/woodsrm/SimTools/pynbody/pynbody/analysis/angmom.py in faceon(h, *_kwargs)
98 the halo of interest), this routine centers the simulation and
99 rotates it so that the disk lies in the x-y plane. This gives
100 a face-on view for SPH images, for instance."""
101
--> 102 sideon(h, calc_faceon_matrix, *_kwargs)

/home/woodsrm/SimTools/pynbody/pynbody/analysis/angmom.py in sideon(h, vec_to_xform, cen_size, disk_size, cen, vcen, verbose)
62 if verbose :
63 print "Finding halo velocity center..."
---> 64 cen = h.star[filt.Sphere(cen_size)]
65 if len(cen)<5 :
66 # fall-back to DM

/home/woodsrm/SimTools/pynbody/pynbody/snapshot.py in getattr(self, name)
268
269 try:
--> 270 return self[family.get_family(name)]
271 except ValueError :
272 pass

/home/woodsrm/SimTools/pynbody/pynbody/snapshot.py in getitem(self, i)
96 return SubSnap(self, i)
97 elif isinstance(i, family.Family) :
---> 98 return FamilySubSnap(self, i)
99 elif isinstance(i, (list, tuple, np.ndarray, filt.Filter)) :
100 return IndexedSubSnap(self, i)

/home/woodsrm/SimTools/pynbody/pynbody/snapshot.py in init(self, base, fam)
930 # Use the slice attributes to find sub array length

931         self._num_particles = self._slice.stop-self._slice.start

--> 932 self._file_units_system = base._file_units_system
933
934 def delitem(self, name) :

/home/woodsrm/SimTools/pynbody/pynbody/snapshot.py in getattr(self, name)
272 pass
273
--> 274 raise AttributeError("%r object has no attribute %r"%(type(self).name, name))
275
276

AttributeError: 'Halo' object has no attribute '_file_units_system'

Original issue: http://code.google.com/p/pynbody/issues/detail?id=36

Original author: [email protected] (January 19, 2011 22:28:38)

Guys: The code is getting complicated!!!!

I just tried to open a simple tipsy file and got a simple tipsy file and got an error that point to a place with so many underscores that I have no idea what is going on:

>>> import pynbody, pynbody.plot
>>> s = pynbody.load('g25271.00880.groupone')
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/shadow/data1/gss/home/src/pynbody/pynbody/init.py", line 24, in load
if c._can_load(filename) : return c(filename,args,*kwargs)
File "/shadow/data1/gss/home/src/pynbody/pynbody/tipsy.py", line 79, in init
self_type[name][n_done:n_done+n_block] = g[:,i]
File "/shadow/data1/gss/home/src/pynbody/pynbody/snapshot.py", line 73, in getitem
raise KeyError(i)
KeyError: 'x'

It's unclear to me why it is having trouble loading "x". Seems pretty straightforward and implemented in tipsy.py.

It looks like snapshot.py is trying to fill in an array that isn't created yet. need to add call the _create_item or _set_item?

But seriously, this code is getting messy. It'd better be doing something awesome useful.

Full disclosure:
a) groupone is a tipsy file of the major halo in a renormalized run.
b) When I just updated, there was some merging going on, it threw me up into a 3 panelled vim window and I no idea what to do, so I just wq'd a bunch times.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=23

Original author: [email protected] (March 21, 2011 19:48:56)

a) It seems to me that the 'time' property should be assigned units that don't have an expansion factor in them like it is now:

In [7]: s.properties['time']
Out[7]: Unit("1.40e+01 s kpc a km**-1")

b) This may be the way we want things, but at the moment, physical_units doesn't touch the properties:

In [5]: s.physical_units()

c) Converting the above units fails:

In [6]: s.properties['time'].in_units('yr')
ERROR: An unexpected error occurred while tokenizing input
The following traceback may be corrupted or invalid
The error message is: ('EOF in multi-line statement', (53, 0))

ERROR: An unexpected error occurred while tokenizing input
The following traceback may be corrupted or invalid
The error message is: ('EOF in multi-line statement', (58, 0))

UnitsException Traceback (most recent call last)

/shadow/data2/gsstinson/g15784/bigstarlog/01024/<ipython console> in <module>()

/shadow/data1/gss/home/src/pynbody/pynbody/units.pyc in in_units(self, _a, *_kw)
182 """Alias for ratio"""
183
--> 184 return self.ratio(_a, *_kw)
185
186 def irrep(self) :

/shadow/data1/gss/home/src/pynbody/pynbody/units.pyc in ratio(self, other, **substitutions)
177 return (self/other).dimensionless_constant(**substitutions)
178 except UnitsException :
--> 179 raise UnitsException, "Not convertible"
180
181 def in_units(self, _a, *_kw) :

UnitsException: Not convertible

In [7]: s.properties['time']
Out[7]: Unit("1.40e+01 s kpc a km**-1")

Original issue: http://code.google.com/p/pynbody/issues/detail?id=43

Original author: [email protected] (March 22, 2011 18:55:08)

It is hard to make meaningful errors come out of the lazy loading __getitem__code as it is now.

At the moment, I want the feedback when one tries s['massform'] to be more informative. If I raise IOError, "Cannot find starlog file" from load_starlog, the IOError gets caught and the less informative KeyError: Unknown array massform is presented to the user.

Would it be possible to switch the order of loading derived arrays and arrays from files?

Rok sort of alluded to this in his recent addition of sound speed elements (which is awesome!), but he has an issue raising an IOError from a derived array. In that case, it seems crucial for the user to be told that the dsp file can't be found.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=44

Original author: [email protected] (January 31, 2011 15:49:13)

What steps will reproduce the problem?

In [3]: s = pynbody.load('g15784.01024')
In [4]: h = s.halos()
In [5]: pynbody.plot.rotation_curve(h[1])

What is the expected output? What do you see instead?

Expected: Pretty plot

Instead:

TypeError Traceback (most recent call last)

/2/home/woodsrm/Data/Gas/g15784.new/<ipython console> in <module>()

/home/woodsrm/SimTools/pynbody/pynbody/plot/profile.py in rotation_curve(s, center, r_units, v_units, disk_height, nbins, bin_spacing)
14
15 if center :
---> 16 angmom.faceon(s)
17
18

/home/woodsrm/SimTools/pynbody/pynbody/analysis/angmom.py in faceon(h)
69 the halo of interest), this routine centres the simulation and
70 rotates it so that the disk lies in the x-y plane. This gives
71 a face-on view for SPH images, for instance."""
72
---> 73 sideon(h, calc_faceon_matrix)

/home/woodsrm/SimTools/pynbody/pynbody/analysis/angmom.py in sideon(h, vec_to_xform, cen_size, disk_size)
57 cen = h.gas[filt.Sphere(disk_size)]
58
---> 59 trans = vec_to_xform(ang_mom_vec(cen))
60
61 top.transform(trans)

/home/woodsrm/SimTools/pynbody/pynbody/analysis/angmom.py in ang_mom_vec(snap)
6
7 The return units are [mass][dist][vel] as per the units of the snapshot."""
----> 8 angmom = (snap['mass'].reshape(len(snap),1)*np.cross(snap['pos'], snap['vel'])).sum(axis=0).view(np.ndarray)
9 return angmom
10

/home/woodsrm/SimTools/pynbody/pynbody/array.py in q(s, _y, *_kw)
532 def _wrap_fn(w) :
533 def q(s, _y, *_kw) :
--> 534 return w(SimArray(s), _y, *_kw)
535 q.name = w.name
536 return q

TypeError: invalid keyword argument

Please use labels and text to provide additional information.

This is just the code from the front page. I have the newest revisions of the code.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=25

Original author: rokroskar (March 09, 2011 13:25:05)

>>> s.computable_keys()
['nefe', 'hetot', 'ofe', 'ne', 'feh', 'mu', 'p', 'mgfe', 'u', 'HeIII', 'sife', 'HII', 'hydrogen', 'vrxy', 'ke', 'smooth', 'j', 'vcxy', 'v2', 'rxy', 'vr', 'r', 'rho']

Note that pynbody is telling me I can get 'mgfe' and 'sife' computed for this snapshot.

>>> s['mgfe']
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/home/itp/roskar/nbody/python/pynbody/pynbody/snapshot.py", line 94, in getitem
raise KeyError, "Unknown array "+i
KeyError: 'Unknown array mgfe'

Of course it doesn't work because my outputs don't have the data for those species.

@TipsySnap.derived_quantity
def mgfe(self) :

in tipsy.py should check whether or not the required arrays are available, and if they are not, don't add the method to the derived_quantity_registry

There should probably be a general mechanism that does this at load time for all derived_quantities - is the best place to do it in snapshot.py:

@classmethod
def derived_quantity(cl,fn):
if not SimSnap._derived_quantity_registry.has_key(cl) :
SimSnap._derived_quantity_registry[cl] = {}
SimSnap._derived_quantity_registry[cl][fn.name]=fn
return fn

I guess the trouble is that you would have to have a dependency list for each derived_quantity

Original issue: http://code.google.com/p/pynbody/issues/detail?id=33

Original author: rokroskar (December 08, 2010 11:30:48)

>>> s = pynbody.load(file)
>>> p = s[0] #doesn't work
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/home/itp/roskar/nbody/python/pynbody/pynbody/snapshot.py", line 65, in getitem
raise TypeError
TypeError
>>> p = s[0,1] # works
>>> p
<SimSnap "file:indexed" len=2>

Original issue: http://code.google.com/p/pynbody/issues/detail?id=8

Original author: [email protected] (March 23, 2011 18:53:35)

Would be useful to have cube selection filters, equivalent to

np.where((f['x']>-xa)(f['x']<xa)(f['y']>-ya)(f['y']<ya)(f['z']>-za)*(f['z']<za))

Original issue: http://code.google.com/p/pynbody/issues/detail?id=46

Original author: [email protected] (March 16, 2011 18:52:03)

It would be nice to be able to control the radius that is plotted in the "rotation curve" routine. Just 2 parameters like "r_min" and "r_max"?

Original issue: http://code.google.com/p/pynbody/issues/detail?id=39

Original author: rokroskar (January 13, 2011 14:42:53)

currently only ascii format is supported for arrays read from auxilliary files. It would be nice to also to support binary format.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=21

Original author: [email protected] (March 15, 2011 21:28:03)

What steps will reproduce the problem?
In [1]: import pynbody
In [2]: import pynbody.analysis.angmom as angmom
In [3]: import pynbody.plot.sph as sph
In [4]: s = pynbody.load('g15784.00168')
In [5]: s.physical_units()
In [6]: h = s.halos()
In [7]: pynbody.plot.rotation_curve(h[1])

What is the expected output? What do you see instead?
expected: plot of rotation curve
instead:

In [23]: pynbody.plot.rotation_curve(h[1])

AttributeError Traceback (most recent call last)

/2/home/woodsrm/Data/Gas/g15784.new/<ipython console> in <module>()

/home/woodsrm/SimTools/pynbody/pynbody/plot/profile.py in rotation_curve(s, center, r_units, v_units, disk_height, nbins, bin_spacing)
20 pro = profile.Profile(s, nbins = nbins, type=bin_spacing)
21
---> 22 r = pro.rbins.in_units(r_units)
23 v = pro['v_circ'].in_units(v_units)
24

AttributeError: Profile instance has no attribute 'rbins'

Original issue: http://code.google.com/p/pynbody/issues/detail?id=38

Original author: [email protected] (December 02, 2010 16:28:39)

numpy.loadtxt is extraordinarily slow. Its use in tipsy.py to read ASCII arrays is incredibly frustrating. We should replace it as soon as possible.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=7

Original author: rokroskar (December 13, 2010 14:05:55)

it would be nice to be able to convert any quantity you have back to the original system units with a single command, i.e.

>>> s['pos'].convert_units('cm') -> converts from whatever units we're currently in to cm
>>> s['pos'].convert_units('system') -> throw it back to the original simulation units

Original issue: http://code.google.com/p/pynbody/issues/detail?id=13

Original author: [email protected] (November 12, 2010 07:11:23)

It would be pretty easy to create keys of the object for all the files that start with the filename. Load the data from the files when needed.
['HI']
['coolontime']
['grp']
['z1.345.amiga_halos']

Can we be intelligent about only loading the data into the particle type that need them?

Original issue: http://code.google.com/p/pynbody/issues/detail?id=1

Original author: rokroskar (November 17, 2010 16:49:12)

I like the family system, but one thing seems to be missing for me - there is no way to ask a SimSnap? instance which families it has defined. It seems to me that the over-arching simulation object should know which families are present, i.e. I want to be able to call s.families() just like I can do s.keys() - or is there some reason why this isn't done?

Original issue: http://code.google.com/p/pynbody/issues/detail?id=2

Original author: rokroskar (December 13, 2010 11:07:29)

accessing the 'r' array calculates the radius from 'pos' and stores it for later use.

>>> s.s['r']
SimArray([ 1.90365715, 1.95675759, 4.49171212, ..., 6.47421208,
10.23414716, 10.68530362], 'kpc')

>>> s['r']
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/Users/rokstar/Nbody/python/pynbody/pynbody/snapshot.py", line 48, in getitem
self._assert_not_family_array(i)
File "/Users/rokstar/Nbody/python/pynbody/pynbody/snapshot.py", line 250, in _assert_not_family_array
raise KeyError, "Array "+name+" is a family-level property"
KeyError: 'Array r is a family-level property'

This should instead generate the 'r' array for the families that don't yet have it, if this is possible.

I can generate 'r' is for the other two families by hand:

>>> s.d['r']
SimArray([ 245.4951134 , 6.18237036, 4.86145924, ..., 372.86583608,
186.47563762, 102.27377937], 'kpc')
>>> s.g['r']
SimArray([ 150.93175958, 67.12237074, 7.0066101 , ..., 19.53358571,
106.83367786, 35.55015996], 'kpc')

now array 'r' is defined for all particle families, but still:

>>> s['r']
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/Users/rokstar/Nbody/python/pynbody/pynbody/snapshot.py", line 48, in getitem
self._assert_not_family_array(i)
File "/Users/rokstar/Nbody/python/pynbody/pynbody/snapshot.py", line 250, in _assert_not_family_array
raise KeyError, "Array "+name+" is a family-level property"
KeyError: 'Array r is a family-level property'

Original issue: http://code.google.com/p/pynbody/issues/detail?id=12

Original author: [email protected] (March 23, 2011 18:53:36)

Would be useful to have cube selection filters, equivalent to

np.where((f['x']>-xa)(f['x']<xa)(f['y']>-ya)(f['y']<ya)(f['z']>-za)*(f['z']<za))

Original issue: http://code.google.com/p/pynbody/issues/detail?id=47

Original author: [email protected] (March 08, 2011 02:53:05)

What steps will reproduce the problem?
In [3]: s = pynbody.load('g15784.00256')

In [4]: h = s.halos()

In [5]: pynbody.plot.rotation_curve(h[1])

What is the expected output? What do you see instead?
Expected: pretty plot

Received:

AttributeError Traceback (most recent call last)

/2/home/woodsrm/Data/Gas/g15784.new/<ipython console> in <module>()

/home/woodsrm/SimTools/pynbody/pynbody/plot/profile.py in rotation_curve(s, center, r_units, v_units, disk_height, nbins, bin_spacing)
20 pro = profile.Profile(s, nbins = nbins, type=bin_spacing)
21
---> 22 r = pro.r.in_units(r_units)
23 v = pro['v_circ'].in_units(v_units)
24

AttributeError: Profile instance has no attribute 'r'

Please use labels and text to provide additional information.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=30

Original author: [email protected] (March 10, 2011 13:02:54)

In agreement with our decision to go with American naming conventions, shouldn't catalogue always actually be catalog?

Wanted to check with the Americans amongst us before making that change.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=34

Original author: rokroskar (December 13, 2010 17:49:31)

I can't figure this one out...

>>> s = pynbody.load('10/12M_hr.01000')
Loading units from /lustre/zbox3b/roskar/isolated_runs/12M_hr/12M_hr.param
>>> s.star['smooth']
SimArray([ NaN, NaN, NaN, ..., NaN, NaN, NaN], 'kpc')
>>> s.star.calculate_smoothing()
>>> s.star['smooth']
SimArray([ 0.00699206, 0.01383823, 0.03779312, ..., 0.04534013,
0.10665303, 0.13607457], 'kpc')

I don't know why the lazy-load isn't working properly, please advise!

Original issue: http://code.google.com/p/pynbody/issues/detail?id=16

Original author: [email protected] (December 13, 2010 15:19:55)

At the moment, Rok's calculate_smoothing will give different results depending on what it's called on. So, for instance

f['gas']['hs']

will calculate a bona-fide gas smoothing length, storing in hs, whereas if someone has first issued

f['hs']

the snapshot-level array already has stored a hybrid all-particle smoothing length.

It's useful, of course, to be able to calculate smoothing lengths for either

(a) one -particle cases
OR (b) all-particle cases

Perhaps we need a naming convention to distinguish this? So

f.gas['family_smooth']

means something different from

f['all_smooth']

As a final related point, the SPH code is using the array name 'smooth' whereas the smoothing length calculation is generating something called 'hs'. This name inconsistency should be fixed. Personally I prefer 'smooth' -- it's less cryptic -- but I don't mind particularly. Obviously the way we resolve this is linked to how we distinguish different modes of smoothing length calculation.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=14

Original author: [email protected] (December 01, 2010 23:20:49)

always does byte swap whether it needs to or not.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=6

Original author: [email protected] (December 10, 2010 11:08:52)

We should implement pynbody.plot.rotation_curve so that our front-page example generates a gorgeous plot every time ;)

Rok?

Original issue: http://code.google.com/p/pynbody/issues/detail?id=10

Original author: rokroskar (December 08, 2010 13:24:28)

I've implemented a simple function called rotateup in snapshot.py:

def rotateup(self, angle):
"""Rotates the snapshot along the current x-axis by 'angle' degrees."""
angle *= np.pi/180
self.transform(np.matrix([[1, 0, 0],
[0, np.cos(angle), -np.sin(angle)],
[0, np.sin(angle), np.cos(angle)]]))

using this on a snapshot 's':
>>> s.rotateup(90)
and plotting it
>>> plot(s['x'],s['y'], ',')

yields what you would expect... but doing

>>> pynbody.plot.sph.image(s.gas)

plots the system with x and y axes inverted... any idea what's going on?

Original issue: http://code.google.com/p/pynbody/issues/detail?id=9

Original author: rokroskar (December 01, 2010 15:49:54)

the file loader should not attempt to do anything with units from the .param file if this is not requested... in particular, if simple (isolated) runs are being analyzed, there is no reason to have a hubble constant specified etc. There should be a flag called "set_units" (or something similar) in the load function that is False by default, but can be turned on by the user.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=5

Original author: [email protected] (January 05, 2011 19:53:47)

I just installed the python.org 2.6.6 pkg followed by the scipy 0.8 Mac OS X pkg and I'm getting a seg fault when I try to run the code on the front page. I tried grabbing the newest numpy 1.5.1 pkg and that didn't seem to help.

What steps will reproduce the problem?

kemosabe:12M stinson$ python
Python 2.6.1 (r261:67515, Jun 24 2010, 21:47:49)
[GCC 4.2.1 (Apple Inc. build 5646)] on darwin
Type "help", "copyright", "credits" or "license" for more information.
>>> import pynbody, pynbody.plot
>>> import pynbody.plot.sph
>>> s = pynbody.load('12M.00350')
Loading units from /Users/stinson/12M/12M.param
>>> pynbody.plot.rotation_curve(s)

Original issue: http://code.google.com/p/pynbody/issues/detail?id=19

Original author: rokroskar (February 11, 2011 15:30:59)

with some outputs I get the following error:

Traceback (most recent call last):
File "/home/itp/roskar/nbody/python/delta_cofm.py", line 13, in <module>
s = pynbody.load(file)
File "/home/itp/roskar/nbody/python/pynbody/pynbody/init.py", line 24, in load
if c._can_load(filename) : return c(filename,args,*kwargs)
File "/home/itp/roskar/nbody/python/pynbody/pynbody/tipsy.py", line 80, in init
g = np.fromstring(f.read(len(st)_n_block_4),'f').byteswap().reshape((n_block,len(st)))
ValueError: total size of new array must be unchanged

There is nothing funny that I can immediately find with the output...

Original issue: http://code.google.com/p/pynbody/issues/detail?id=26

Original author: [email protected] (March 23, 2011 19:21:00)

This is killing me:

Package scipy-0.7.0-5.fc12.x86_64 already installed and latest version

scipy-0.9 has a whole new implementation of 2d interpolation that seems to do things the much more sensible way.

This may be me crying in the wilderness. The interpolation suggests using map_coordinates:

http://www.scipy.org/Cookbook/Interpolation#head-6a59d9e4e9645b568f7877809306bff1130fa34b

I will see if that works any better in the morning.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=48

Original author: [email protected] (January 06, 2011 08:45:15)

Subviews are a critical component of the pynbody model, e.g. one can load halo catalogues and request just the particles in halo 1.

f = pynbody.load(...)
h_1 = f.halos()[1]

When accessing a non-existent array from one of these subviews, it is possible that the non-existent array name will be matched to a 'derivable'. For the purposes of this example, let's take the example of 'r', which in derived.py is defined to be the radial distance from the origin.

If we access f['r'], then all the radial positions are calculated. However, if we had first accessed h_1['r'], the array would be created but only calculated for the particles within halo 1. The rest of the array is filled with zeros. Then the derived_quantity routine is not called again, because the array already exists.

Therefore, depending on previous operations, f['r'] may return different/wrong information. There is a mechanism to avoid this, by calling f.derive('r'), but this seems rather wasteful.

What is the correct solution? Various possibilities present themselves

1. Derived arrays are NOT stored, unless specifically requested. So asking for h_1['r'] returns the correct information but recalculates it each time it is requested unless one calls h_1['r'].store() or some such. This rather negates the idea of the derived arrays, which was to have a transparent way of getting information that may or may not have yet been calculated (since we need to be aware of the state of the returned array.)

Or 2. We have some notion of arrays which are local to a subview, so we can store h_1['r'] locally to h_1 without ever creating the 'r' array in f. The trouble here is that the h_1['r'] array can now mask a legitimately created 'r' array in f.

Or 3. We have some notion of masks on arrays, so that 'uncalculated' values can be flagged. But this introduces annoying conceptual and computational overheads.

Original issue: http://code.google.com/p/pynbody/issues/detail?id=20

Original author: [email protected] (January 05, 2011 18:18:51)

What steps will reproduce the problem?

import pynbody
s = pynbody.load('filename')
print s.properties
> {'a': 0.33356858444748205, 'h': 0.72987934710075975, 'omegaL0': '0.76', 'z': 1.9978842331821651, 'omegaM0': '0.24'}

What is the expected output? What do you see instead?

Expected: {'a': 0.33356858444748205, 'h': 0.72987934710075975, 'omegaL0': 0.76, 'z': 1.9978842331821651, 'omegaM0': 0.24}

Please use labels and text to provide additional information.

In tipsy.py - no "float()" call on omegaM0 or omegaL0.

    try :
        s = line.split()
        if s[0]=="dMsolUnit" :
        munit_st = s[2]+" Msol"
        munit = float(s[2])
        elif s[0]=="dKpcUnit" :
        dunit_st = s[2]+" kpc a"
        dunit = float(s[2])
        elif s[0]=="dHubble0" :
        hub = float(s[2])
        elif s[0]=="dOmega0" :
        om_m0 = s[2]
        elif s[0]=="dLambda" :
        om_lam0 = s[2]

Original issue: http://code.google.com/p/pynbody/issues/detail?id=18

Original author: [email protected] (March 14, 2011 16:31:30)

What steps will reproduce the problem?
In [1]: import pynbody, pynbody.plot.sph as sph, pynbody.analysis.angmom as angmom
In [5]: s = pynbody.load('g15784.00168')
In [8]: sph.image(s.gas,qty='PdV',width=200)

What is the expected output? What do you see instead?
Expected: Pretty plot.

Instead:

ValueError Traceback (most recent call last)

/2/home/woodsrm/Data/Gas/g15784.new/<ipython console> in <module>()

/home/woodsrm/SimTools/pynbody/pynbody/plot/sph.py in image(sim, qty, width, resolution, units, log, vmin, vmax, av_z)
75
76 if log :
---> 77 im[np.where(im==0)] = abs(im[np.where(im!=0)]).min()
78 im = np.log10(im)
79

ValueError: zero-size array to ufunc.reduce without identity

Please use labels and text to provide additional information.

I have run this before and gotten an image out. The array looks like:

In [6]: max(s.g['PdV'])
Out[6]: 3816.9532388748999
In [7]: min(s.g['PdV'])
Out[7]: -4321.5720493875997

I negatives didn't used to cause a problem, so I don't think they are now?

Original issue: http://code.google.com/p/pynbody/issues/detail?id=35