OH! And what's the radius of the halo, and how many cells does it get if you make a sphere from the HOPGroup? THAT is the key thing, actually... On Sun, Jan 25, 2009 at 3:20 PM, Matthew Turk <matthewturk@gmail.com> wrote:

Okay, so this could in fact be a problem with the Profile.py module.

Here's the info I'd need -- what are the radii it's trying to calculate the profile over, and what's the number of bins it's trying to do? How many bins are you generating?

-Matt

On Sun, Jan 25, 2009 at 3:15 PM, Stephen Skory <stephenskory@yahoo.com> wrote:

...

Hi,

I'm having some trouble with the HaloProfiler. I'm running with the trunk r1147 sample_halo_profile.par and runHaloProfiler.py. Here's the error I'm getting:

File "runHaloProfiler.py", line 8, in <module> q.makeProfiles() File "/share/home/00649/tg457850/yt/lib/python2.6/site-packages/yt-1.5dev.dev_r1147-py2.6-linux-x86_64.egg/yt/extensions/HaloProfiler.py", line 107, in makeProfiles virial = self._CalculateVirialQuantities(profile) File "/share/home/00649/tg457850/yt/lib/python2.6/site-packages/yt-1.5dev.dev_r1147-py2.6-linux-x86_64.egg/yt/extensions/HaloProfiler.py", line 276, in _CalculateVirialQuantities if (overDensity[1] <= self.haloProfilerParameters['VirialOverdensity']): IndexError: list index out of range

For some of the haloes, profile['CellVolume'] and profile['TotalMassMsun'] is being set to all zeros. If I put this 'if' below in HaloProfiler.py before the call to CalculateVirialQuantities(), say line 104, I don't get the error above and the script runs to completion, making profiles and images.

if profile['CellVolume'][0] == 0.0 : continue

profile['CellVolume'] being zero gives nans in _AddActualOverdensity() in HaloProfiler.py.

If I run this in a unaltered enzo dataset (with no .yt, HopAnalysis.out etc... files), I get this problem for every dataset I've tried. I've tried two runs of 64^3 root grid with 3 levels cosmology, one is 1 mpc/h and the other 64 mpc/h.

I've done some sleuthing, but I don't undertstand enough of how Profiles.py to make sense of it. 'args' in _lazy_add_fields() is returning empty from _get_bins(), and _get_bins() is returning empty because self._get_field() is returning nothing.

Those of you who understand this better, how can I help you?

Thanks!

_______________________________________________________ sskory@physics.ucsd.edu o__ Stephen Skory http://physics.ucsd.edu/~sskory/ _.>/ _Graduate Student ________________________________(_)_\(_)_______________ _______________________________________________ Yt-dev mailing list Yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org

I think I see what the problem is here. The smallest bin for the profile is 0.021 and the maximum radius is 0.029, so even if 64 bins are requested, the profiler is only going to have one bin or bins full of zeros. From the orginal error message, it looks like it only had one bin. Since virial quantities are calculated by interpolating between bins where the overdensity crosses the critical overdensity, it needs at least 2 bins to do anything. I think the solution is to put a check in to see if the profile only has one bin, and skip further operations on that halo. As a side note, the first bin for CellVolume will always be zero, since the radius of the bin is zero. Maybe that should be changed. Any thoughts on that? I can commit a fix to this that works as mentioned above. Any objections to this? Thanks for finding that! Britton On Sun, Jan 25, 2009 at 4:27 PM, Stephen Skory <stephenskory@yahoo.com>wrote:

----- Original Message ----

...

From: Matthew Turk <matthewturk@gmail.com> To: yt-dev@lists.spacepope.org Sent: Sunday, January 25, 2009 3:22:09 PM Subject: Re: [Yt-dev] HaloProfiler Problems

OH! And what's the radius of the halo, and how many cells does it get if you make a sphere from the HOPGroup? THAT is the key thing, actually...

When it calls this on the first problem halo:

profile = lagos.BinnedProfile1D(sphere,self.haloProfilerParameters['n_bins'],"RadiusMpc", (4*self.pf.h.get_smallest_dx() * self.pf ['mpc']),halo['r_max'], log_space=True, lazy_reader=True)

self.haloProfilerParameters['n_bins'] = 64 (4*self.pf.h.get_smallest_dx() * self.pf['mpc']) = 0.021143 halo['r_max'] = 0.029493

_______________________________________________________ sskory@physics.ucsd.edu o__ Stephen Skory http://physics.ucsd.edu/~sskory/ <http://physics.ucsd.edu/%7Esskory/> _.>/ _Graduate Student ________________________________(_)_\(_)_______________ _______________________________________________ Yt-dev mailing list Yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org

Thanks for finding that!

You're welcome. I guess you're just used to running your code on actual high-resolution science-quality datasets, while I'm just playing around with crummy rough test cases. _______________________________________________________ sskory@physics.ucsd.edu o__ Stephen Skory http://physics.ucsd.edu/~sskory/ _.>/ _Graduate Student ________________________________(_)_\(_)_______________

I just tried r1148 and I got the same problem as before. Sorry for the slow response. _______________________________________________________ sskory@physics.ucsd.edu o__ Stephen Skory http://physics.ucsd.edu/~sskory/ _.>/ _Graduate Student ________________________________(_)_\(_)_______________

It's below. The only thing I've done differently from the sample provided is comment out HopOutputFile (so it uses the default HopAnalysis.out), and I changed the value of ProjectionWidth, but I don't think that has any effect on the issue at hand. # Sample HaloProfiler parameter file. ProfileOutputDir = radial_profiles ProjectionOutputDir = projections # Name of hop output file #HopOutputFile = s-1143.txt # Name of output file with virial quantities VirialQuantitiesOutputFile = VirialQuantities.out # Minimum virial mass in solar masses VirialMassCutoff = 1e7 # Critical overdensity for virial cutoff. VirialOverdensity = 200 n_bins = 64 # 1D radial profiles # Syntax is Profile[Field] = Weight_Field,Accumulate (True or False) Profile[CellVolume] = None,True Profile[TotalMassMsun] = None,True # Projections of individual halos # Syntax is Projection[Field] = Weight_Field Projection[Density] = None Projection[Temperature] = Density ProjectionWidth = 1.0 # Mpc # Max level at which halos will be projected. # This determines the size of the fixed resolution buffer. # Can be set to any level or use the word 'max' to do the # max level in the simulation. ProjectAtLevel = max _______________________________________________________ sskory@physics.ucsd.edu o__ Stephen Skory http://physics.ucsd.edu/~sskory/ _.>/ _Graduate Student ________________________________(_)_\(_)_______________

Yes, we know. Britton's just committed a fix... It's been tested on "RD0005-mine." On Mon, Jan 26, 2009 at 9:50 AM, Stephen Skory <stephenskory@yahoo.com> wrote:

All,

+ virial = {} + if (len(profile[overdensity_ field]) < 2): + for field in fields: + virial[field] = 0.0 + return virial

This isn't going to fix it. 'print profile[overdensity_field]' gives this for the bad halo:

[ NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN]

Again, sorry for not following up more quickly on this.

_______________________________________________________ sskory@physics.ucsd.edu o__ Stephen Skory http://physics.ucsd.edu/~sskory/ _.>/ _Graduate Student ________________________________(_)_\(_)_______________ _______________________________________________ Yt-dev mailing list Yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org

Let me know if this comes up again.

I will, but it appears fixed on my end too. Thanks Britton and Matt! _______________________________________________________ sskory@physics.ucsd.edu o__ Stephen Skory http://physics.ucsd.edu/~sskory/ _.>/ _Graduate Student ________________________________(_)_\(_)_______________