[yt-users] Hop_particles

Stephen Skory stephenskory at yahoo.com
Mon Nov 30 11:16:54 PST 2009


Shankar,

> I mean, instead of using...
> 
> p.modify["hop_particles"](halos, p_size=0.05, max_number=4000, min_size=0, 
> alpha=0.6)
> 
> Can I use...
> 
> p.modify["hop_particles"]("HopAnalysis.out", p_size=0.05, max_number=4000, 
> min_size=0, alpha=0.6)
> 
> 
> I tried it but did not work. I had to generate the halos again.

Python is a programming language. Therefore, you cannot just put strings in place of other kinds of items and hope that it will work. In your code above, 'halos' is a HaloList object, which as a kind of Python construct that is defined in yt for handling a list of halo data objects. You can see where it's defined (look for "class HaloList(object):") in yt/lagos/HaloFinding.py.

To kind of answer your real question, there is a way to plot halo circles without running HOP again. I have previously helped Carolyn Peruta with this (but I didn't send it to this list). Below is what I sent her. However, this doesn't help you because you want to plot hop_particles from a previous run of HOP. If you look at your HopAnalysis.out file, you'll see that there is no way to recover the identities of particles in a halo from that simple text file. Remember that a halo isn't usually spherical, so a center of mass and radius is not enough to get the true particle membership.

The way you'd have to do it is to save the particles to a HDF5 file and read in again when plotting. You can write the particle data to HDF5 files using "halos.write_particle_list('prefix')" (which saves to a file named prefix.h5, or files named prefixNNNN.h5, NNNN depending on how many tasks you're using). There is a companion function "halos.write_particle_list_txt('prefix')" which writes a text file which lists in which .h5 file the particle data is for each halo (this is only required when running with more than one task).

As a quick sketch, here is how one might recover the positions of particles from all the haloes from one .h5 file:

import h5py
fn = h5py.File('halos.h5')
old_haloes_parts = []
for halo_group in fn:
    xpos = fn[halo_group]["particle_position_x"][:]
    ypos = fn[halo_group]["particle_position_y"][:]
    zpos = fn[halo_group]["particle_position_z"][:]
    old_haloes_parts.append((xpos,ypos,zpos))

then you could feed old_haloes_parts along with a old_haloes (as constructed below), to a modified "hop_particles" callback and plot the particles. Please study the modified "hop_circles" below to see how to change "hop_particles". I'm sorry I don't have time to do this for you, but I think I've given you a fair start. Good luck!

(p.s. Adding a built-in function to re-read in haloes is something I'm considering doing. But it would be different than how it's done below (to be more general), which will require a bit more thought and time on my part.)

----

This first stuff goes in the same file that you drive yt with. It can go before or after
the yt stuff because it's independent of yt.

What I'm doing is just making a list where each entry contains the essential bits
for each halo. This list is then passed to the slightly modified OldHopCircleCallback
that expects a list of these things rather than a list of Halo objects.

You'll have to fix the columns for the values of the halo quantities below,
below the 'while line:' line. Hopefully I am not patronizing you by reminding you
that Python is zero-ordered, so the first column is 0.

fp = open("HopAnalysis.out", "r")
old_haloes = []
line = fp.readline()
# You need one of these lines to start, and then put one of these for each line
# in the text that is a header line. So if there is one header line that gives
# the names of the columns, you shouldn't have to do anything here.
line = fp.readline()
while line:
    line = line.split() # splits the line by spaces
    size = int(line[0]) # Change this to the column for the halo size
    xp = float(line[1]) # column for x position
    yp = float(line[2]) # etc...
    zp = float(line[3])
    mass = float(line[4])
    radius = float(line[5])
    item = [size, xp, yp, zp, radius, mass]
    old_haloes.append(item) # adds the halo to the list
    line = fp.readline() # reads a new line.

fp.close()

(.... do your yt stuff ... )

pc = PlotCollection(pf)
pc.add_projection("Density",0)
pc.plots[-1].add_callback(OldHopCircleCallback(old_haloes))

pc.save("Circles")

------------------------------------

This is added to yt/raven/Callbacks.py.

class OldHopCircleCallback(PlotCallback):
    _type_name = "old_hop_circles"
    def __init__(self, hop_output, max_number=None,
                 annotate=False, min_size=20, max_size=10000000,
                 font_size=8, print_halo_size=False,
                 print_halo_mass=False, width=None):
        self.hop_output = hop_output
        self.max_number = max_number
        self.annotate = annotate
        self.min_size = min_size
        self.max_size = max_size
        self.font_size = font_size
        self.print_halo_size = print_halo_size
        self.print_halo_mass = print_halo_mass
        self.width = width

    def __call__(self, plot):
        from matplotlib.patches import Circle
        x0, x1 = plot.xlim
        y0, y1 = plot.ylim
        l, b, width, height = _get_bounds(plot._axes.bbox)
        xi = lagos.x_dict[plot.data.axis]
        yi = lagos.y_dict[plot.data.axis]
        dx = plot.image._A.shape[0] / (x1-x0)
        dy = plot.image._A.shape[1] / (y1-y0)
        for id, halo in enumerate(self.hop_output[:self.max_number]):
            size = halo[0]
            if size < self.min_size or size > self.max_size: continue
            com = na.array([halo[1], halo[2], halo[3]])
            radius = halo[4]
            tot_m = halo[5]
            if self.width is not None and \
                na.abs(com - 
                       plot.data.center)[plot.data.axis] > \
                   self.width:
                continue
            radius = radius * dx
            center = com
            center_x = (center[xi] - x0)*dx
            center_y = (center[yi] - y0)*dy
            cir = Circle((center_x, center_y), radius, fill=False)
            plot._axes.add_patch(cir)
            if self.annotate:
                if self.print_halo_size:
                    plot._axes.text(center_x, center_y, "%s" % size,
                    fontsize=self.font_size)
                elif self.print_halo_mass:
                    plot._axes.text(center_x, center_y, "%s" % tot_m,
                    fontsize=self.font_size)
                else:
                    plot._axes.text(center_x, center_y, "%s" % id,
                    fontsize=self.font_size)

 _______________________________________________________
sskory at physics.ucsd.edu           o__  Stephen Skory
http://physics.ucsd.edu/~sskory/ _.>/ _Graduate Student
________________________________(_)_\(_)_______________



----- Original Message ----
> From: "Agarwal, Shankar" <sagarwal at ku.edu>
> To: Discussion of the yt analysis package <yt-users at lists.spacepope.org>
> Sent: Sat, November 28, 2009 11:34:24 AM
> Subject: Re: [yt-users] Hop_particles
> 
> Thanks Stephen, that worked. Also, generating halos takes some time. I don't 
> want to run it again and again. Is there a way to directly use HopAnalysis.out 
> to plot hop_circles and hop_particles ? 
> 

> 
> 
> shankar
> 
> 
> -----Original Message-----
> From: yt-users-bounces at lists.spacepope.org on behalf of Stephen Skory
> Sent: Sat 11/28/2009 12:09 PM
> To: Discussion of the yt analysis package
> Subject: Re: [yt-users] Hop_particles
> 
> Shankar,
> 
> 
> >     p.modify["hop_particles"](halos, p_size=1.0, max_number=None, min_size=0, 
> > alpha=0.20000000000000001)
> 
> try setting 'max_number' to something other than 'None', like the number (or 
> greater) of haloes for which you want the particles plotted. With that, it works 
> for me.
> 
> _______________________________________________________
> sskory at physics.ucsd.edu           o__  Stephen Skory
> http://physics.ucsd.edu/~sskory/ _.>/ _Graduate Student
> ________________________________(_)_\(_)_______________
> _______________________________________________
> yt-users mailing list
> yt-users at lists.spacepope.org
> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org




More information about the yt-users mailing list