[yt-users] clump finder

Matthew Turk matthewturk at gmail.com
Wed May 19 15:32:17 PDT 2010


Hi JC,

As a quick note, all this inline text output from IPython is very
difficult to parse over email.  I highly recommend that you script
this task and then run with the --paste command, as mentioned here:

http://yt.enzotools.org/doc/askingforhelp.html

and more specifically, here:

http://yt.enzotools.org/doc/advanced/debugdrive.html#error-reporting

That way, some of the rest of us can try reproducing it, as well.  :)

Anyway, I've looked at the issues, and it looks to me like the first
is because of a mismatch in what the clump callback wants and what the
clump *finder* (as opposed to the contour finder) hands off.  I'll
work with Britton Smith on making that work.

The second error is because the clump finder (again, not the contour
finder) expects to be handed something with a meaningful "Temperature"
field.  This is true for, say, Cosmology simulations and a few types
of star formation simulations.  The Clump finder is a very powerful
tool that also assumes you have relatively rich data; however, you can
mock up a small portion of it by using extract_connected_sets:

http://yt.enzotools.org/doc/modules/amrdata.html#yt.lagos.AMR3DData.extract_connected_sets

and then on each returned object, call the quantity IsBound.  I will
work on writing up a very simple recipe for this, to be included with
the next build of the docs, but Britton Smith or David Collins can
probably chime in with some additional words of wisdom.

As to your third question, IsBound takes as an argument whether or not
to consider the thermal energy.  I believe by default it *does* take
into account thermal energy inside the clump finder, but I defer to
Britton on that.

I hope that helps!

-Matt

On Wed, May 19, 2010 at 3:05 PM, Jean-Claude Passy <jcpassy at gmail.com> wrote:
> Dear all,
>
> I want ultimately plot bound/unbound quantities such as mass, angular
> momentum...
> Britton suggested I used the clump finder. I have read the paragraph in the
> Cookbook but did not manage to make it work. See below:
>
>
> ##########################################################################################
> In [1]: from yt.mods import *
> In [2]: import yt.raven as R
> In [3]: fn = "DD0000/CommonEnvelope0000"
> In [4]: field = "Density"
> In [5]: step = 10.0
> In [6]: pf = load(fn)
> In [7]: whole_box = pf.h.all_data() yt         INFO       2010-05-19
> 13:47:49,272 Getting the binary hierarchy yt         INFO       2010-05-19
> 13:47:49,275 Finished with binary hierarchy reading
> In [8]: total_mass = whole_box.quantities['TotalQuantity']('CellMassMsun')
> In [9]: data_source = pf.h.all_data()
> In [10]: master_clump = Clump(data_source, None, field) yt
> INFO       2010-05-19 13:47:49,371 Getting field Density from 1
> In [11]: c_min = 10**na.floor(na.log10(data_source[field]).min()  )
> In [12]: c_max = 10**na.floor(na.log10(data_source[field]).max()+1)
> In [13]: step = 10
> In [14]: find_clumps(master_clump, c_min, c_max, step) Finding clumps: min:
> 1.000000e-05, max: 1.000000e+03, step: 10.000000 First pass100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 Calculating joins 100%
> |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| Time:
> 00:00:00 yt         INFO       2010-05-19 13:47:49,908 Coalescing 1 joins
> Joining 100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 yt         INFO       2010-05-19 13:47:49,910 Getting field
> tempContours from 1 yt         INFO       2010-05-19 13:47:49,916 Getting
> field tempContours from 1 yt         INFO       2010-05-19 13:47:49,920
> Identified 1 contours between 1.00000e-05 and 1.05005e+02 yt
> INFO       2010-05-19 13:47:49,931 Getting field Density from 1 Finding
> clumps: min: 1.000000e-04, max: 1.000000e+03, step: 10.000000 Wiping out
> existing children clumps. First pass100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 Calculating joins 100%
> |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| Time:
> 00:00:00 yt         INFO       2010-05-19 13:47:49,944 Coalescing 1 joins
> Joining 100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 yt         INFO       2010-05-19 13:47:49,945 Getting field
> tempContours from 1 yt         INFO       2010-05-19 13:47:49,950 Getting
> field tempContours from 1 yt         INFO       2010-05-19 13:47:49,954
> Identified 1 contours between 1.00000e-04 and 1.05005e+02 yt
> INFO       2010-05-19 13:47:49,957 Getting field Density from 1 Finding
> clumps: min: 1.000000e-03, max: 1.000000e+03, step: 10.000000 Wiping out
> existing children clumps. First pass100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 Calculating joins 100%
> |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| Time:
> 00:00:00 yt         INFO       2010-05-19 13:47:49,967 Coalescing 1 joins
> Joining 100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 yt         INFO       2010-05-19 13:47:49,968 Getting field
> tempContours from 1 yt         INFO       2010-05-19 13:47:49,973 Getting
> field tempContours from 1 yt         INFO       2010-05-19 13:47:49,977
> Identified 1 contours between 1.00000e-03 and 1.05005e+02 yt
> INFO       2010-05-19 13:47:49,980 Getting field Density from 1 Finding
> clumps: min: 1.000000e-02, max: 1.000000e+03, step: 10.000000 Wiping out
> existing children clumps. First pass100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 Calculating joins 100%
> |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| Time:
> 00:00:00 yt         INFO       2010-05-19 13:47:49,990 Coalescing 1 joins
> Joining 100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 yt         INFO       2010-05-19 13:47:49,991 Getting field
> tempContours from 1 yt         INFO       2010-05-19 13:47:49,996 Getting
> field tempContours from 1 yt         INFO       2010-05-19 13:47:50,000
> Identified 1 contours between 1.00000e-02 and 1.05005e+02 yt
> INFO       2010-05-19 13:47:50,003 Getting field Density from 1 Finding
> clumps: min: 1.000000e-01, max: 1.000000e+03, step: 10.000000 Wiping out
> existing children clumps. First pass100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 Calculating joins 100%
> |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| Time:
> 00:00:00 yt         INFO       2010-05-19 13:47:50,013 Coalescing 1 joins
> Joining 100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 yt         INFO       2010-05-19 13:47:50,014 Getting field
> tempContours from 1 yt         INFO       2010-05-19 13:47:50,019 Getting
> field tempContours from 1 yt         INFO       2010-05-19 13:47:50,023
> Identified 1 contours between 1.00000e-01 and 1.05005e+02 yt
> INFO       2010-05-19 13:47:50,026 Getting field Density from 1 Finding
> clumps: min: 1.000000e+00, max: 1.000000e+03, step: 10.000000 Wiping out
> existing children clumps. First pass100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 Calculating joins 100%
> |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| Time:
> 00:00:00 yt         INFO       2010-05-19 13:47:50,034 Coalescing 1 joins
> Joining 100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 yt         INFO       2010-05-19 13:47:50,036 Getting field
> tempContours from 1 yt         INFO       2010-05-19 13:47:50,040 Getting
> field tempContours from 1 yt         INFO       2010-05-19 13:47:50,045
> Identified 1 contours between 1.00000e+00 and 1.05005e+02 yt
> INFO       2010-05-19 13:47:50,047 Getting field Density from 1 Finding
> clumps: min: 1.000000e+01, max: 1.000000e+03, step: 10.000000 Wiping out
> existing children clumps. First pass100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 Calculating joins 100%
> |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| Time:
> 00:00:00 yt         INFO       2010-05-19 13:47:50,055 Coalescing 1 joins
> Joining 100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 yt         INFO       2010-05-19 13:47:50,056 Getting field
> tempContours from 1 yt         INFO       2010-05-19 13:47:50,061 Getting
> field tempContours from 1 yt         INFO       2010-05-19 13:47:50,065
> Identified 1 contours between 1.00000e+01 and 1.05005e+02 yt
> INFO       2010-05-19 13:47:50,067 Getting field Density from 1 Finding
> clumps: min: 1.000000e+02, max: 1.000000e+03, step: 10.000000 Wiping out
> existing children clumps. First pass100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 Calculating joins 100%
> |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| Time:
> 00:00:00 yt         INFO       2010-05-19 13:47:50,073 Coalescing 0 joins
> Joining 100%
> ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
> Time: 00:00:00 yt         INFO       2010-05-19 13:47:50,075 Getting field
> tempContours from 1 yt         INFO       2010-05-19 13:47:50,080 Getting
> field tempContours from 1 yt         INFO       2010-05-19 13:47:50,084
> Identified 0 contours between 1.00000e+02 and 1.05005e+02
> In [15]: p = R.PlotCollection(pf) yt.lagos   INFO       2010-05-19
> 13:47:50,918 Max Value is 1.05005e+02 at 0.4843750000000000
> 0.4843750000000000 0.4843750000000000 in grid EnzoGrid_0001 at level 0 (15,
> 15, 15) yt         INFO       2010-05-19 13:47:50,919 Created plot
> collection with default plot-center = [0.484375, 0.484375, 0.484375]
> In [16]: p.add_slice("Density",2,center=[0.5, 0.5, 0.5]) yt
> INFO       2010-05-19 13:47:51,479 Added slice of Density at z = 0.5 with
> 'center' = [0.5, 0.5, 0.5] Out[16]: <yt.raven.PlotTypes.SlicePlot object at
> 0x42f1c10>
> In [17]: p.plots[-1].modify["clumps"](master_clump) Out[17]: 0
> In [18]: p.save('testclumps', format = 'png') 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', (10, 0))
> ---------------------------------------------------------------------------
> TypeError                                 Traceback (most recent call last)
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/scripts/iyt in <module>() ----> 1
>       2       3       4       5
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/raven/PlotCollection.pyc in
> save(self, basename, format, override, force_save)      81
> fn.append(plot.save_image(basename, \      82
> format=format, submit=self._run_id, ---> 83
> override=override, force_save=force_save))      84             if
> self.submit:      85                 im = plot.im.copy()
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/raven/PlotTypes.pyc in
> save_image(self, prefix, format, submit, override, force_save)
> 116         *override* will force no filename generation beyond the prefix.
>     117         """ --> 118         self._redraw_image()     119         if
> not override:     120             self._generate_prefix(prefix)
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/raven/PlotTypes.pyc in
> _redraw_image(self, *args)     404
> self.norm.autoscale(na.array((newmin,newmax)))     405
> self._reset_image_parameters() --> 406         self._run_callbacks()     407
>     408     def _reset_image_parameters(self):
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/raven/PlotTypes.pyc in
> _run_callbacks(self)     259         self._axes.texts = []     260
> for cb in self._callbacks: --> 261             cb(self)     262     263
> def set_label(self, label):
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/raven/Callbacks.pyc in
> __call__(self, plot)     451         nx, ny = plot.image._A.shape
> 452         buff = na.zeros((nx,ny),dtype='float64') --> 453         for
> i,clump in enumerate(reversed(self.clumps)):     454
> mylog.debug("Pixelizing contour %s", i)     455
> TypeError: object of type 'Clump' has no len()
> ##########################################################################################
>
>
> If I also try to reproduce the example in the cookbook:
>
> ##########################################################################################
> In [19]: f = open('%s_clump_hierarchy.txt' % pf,'w')
> In [20]: write_clump_hierarchy(master_clump,0,f) yt         INFO
> 2010-05-19 14:00:24,974 Getting field CellMassMsun from 1 yt
> INFO       2010-05-19 14:00:24,974 Getting field CellVolume from 1
> yt         INFO       2010-05-19 14:00:24,974 Getting field dx from 1
> yt         INFO       2010-05-19 14:00:24,978 Getting field dy from 1
> yt         INFO       2010-05-19 14:00:24,981 Getting field dz from 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', (37, 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', (111, 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', (111, 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', (133, 0))
> ---------------------------------------------------------------------------
> NeedsDataField                            Traceback (most recent call last)
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/scripts/iyt in <module>() ----> 1
>       2       3       4       5
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/Clump.pyc in
> write_clump_hierarchy(clump, level, f_ptr)     196         f_ptr.write("\t")
>     197     f_ptr.write("Clump at level %d:\n" % level) --> 198
> clump.write_info(level,f_ptr)     199     f_ptr.write("\n")     200
> f_ptr.flush()
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/Clump.pyc in
> write_info(self, level, f_ptr)     102                 value =
> item['quantity']()     103             else: --> 104                 value =
> eval(item['quantity'])     105             output = eval(item['format'])
> 106             f_ptr.write("%s%s" % ('\t'*level,output))
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/Clump.pyc in <module>()
> ----> 1       2       3       4       5
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/DerivedQuantities.pyc in
> __call__(self, *args, **kwargs)      70
> self._data_source.pf.h.io)      71         if lazy_reader and not
> self.force_unlazy: ---> 72             return self._call_func_lazy(args,
> kwargs)      73         else:      74             return
> self._call_func_unlazy(args, kwargs)
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/DerivedQuantities.pyc in
> _call_func_lazy(self, args, kwargs)      77         self.retvals = [ [] for
> i in range(self.n_ret)]      78         for gi,g in
> enumerate(self._get_grids()): ---> 79             rv =
> self.func(GridChildMaskWrapper(g, self._data_source), *args, **kwargs)
> 80             for i in range(self.n_ret): self.retvals[i].append(rv[i])
>      81             g.clear_data()
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/DerivedQuantities.pyc in
> _WeightedAverageQuantity(data, field, weight)     151     :param weight: The
> field to weight by     152     """ --> 153     num = (data[field] *
> data[weight]).sum()     154     den = data[weight].sum()     155     return
> num, den
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/DerivedQuantities.pyc in
> __getitem__(self, item)      41         return getattr(self.grid, attr)
> 42     def __getitem__(self, item): ---> 43         return
> self.data_source._get_data_from_grid(self.grid, item)      44      45 class
> DerivedQuantity(ParallelAnalysisInterface):
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/BaseDataTypes.pyc in
> save_state(self, grid, field)      41         old_keys = grid.data.keys()
>      42         grid.field_parameters = self.field_parameters --->
> 43         tr = func(self, grid, field)      44
> grid.field_parameters = old_params      45         grid.data = dict( [(k,
> grid.data[k]) for k in old_keys] )
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/BaseDataTypes.pyc in
> _get_data_from_grid(self, grid, field)    1670         else:
> 1671             pointI = self._get_point_indices(grid) -> 1672
> if grid[field].size == 1: # dx, dy, dz, cellvolume    1673                 t
> = grid[field] * na.ones(grid.ActiveDimensions, dtype='float64')
> 1674                 return t[pointI].ravel()
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/BaseGridType.pyc in
> __getitem__(self, key)     133         """     134         if not
> self.data.has_key(key): --> 135             self.get_data(key)
> 136         return self.data[key]     137
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/BaseGridType.pyc in
> get_data(self, field)     176                     else: raise
> 177             else: --> 178                 self._generate_field(field)
>     179         return self.data[field]     180
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/BaseGridType.pyc in
> _generate_field(self, field)     121                 self[field] =
> temp_array[sl]     122             else: --> 123                 self[field]
> = self.pf.field_info[field](self)     124         else: # Can't find the
> field, try as it might     125             raise exceptions.KeyError, field
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/FieldInfoContainer.pyc in
> __call__(self, data)     316         ii = self.check_available(data)
> 317         original_fields = data.keys() # Copy --> 318         dd =
> self._function(self, data)     319         dd *=
> self._convert_function(data)     320         for field_name in data.keys():
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/UniversalFields.pyc in
> _JeansMassMsun(field, data)     769     770     return (MJ_constant * -->
> 771             ((data["Temperature"]/data["MeanMolecularWeight"])**(1.5)) *
>     772             (data["Density"]**(-0.5)))     773
> add_field("JeansMassMsun",function=_JeansMassMsun,
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/BaseGridType.pyc in
> __getitem__(self, key)     133         """     134         if not
> self.data.has_key(key): --> 135             self.get_data(key)
> 136         return self.data[key]     137
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/BaseGridType.pyc in
> get_data(self, field)     176                     else: raise
> 177             else: --> 178                 self._generate_field(field)
>     179         return self.data[field]     180
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/BaseGridType.pyc in
> _generate_field(self, field)     111             # First we check the
> validator
>     112             try: --> 113
> self.pf.field_info[field].check_available(self)     114             except
> NeedsGridType, ngt_exception:     115                 # This is only going
> to be raised if n_gz > 0
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/FieldInfoContainer.pyc in
> check_available(self, data)     283         """     284         for
> validator in self.validators: --> 285             validator(data)
> 286         # If we don't get an exception, we're good to go
>     287         return True
> /rpod2/jcpassy/yt-x86_64/src/yt-trunk-svn/yt/lagos/FieldInfoContainer.pyc in
> __call__(self, data)     381                 doesnt_have.append(f)
> 382         if len(doesnt_have) > 0: --> 383             raise
> NeedsDataField(doesnt_have)     384         return True     385
> NeedsDataField: (['Temperature'])
> ##########################################################################################
>
>
> Also, if I correctly understand the Cookbook and the reference, the tree of
> clumps in contained in master_clump. It contains only the bound clumps so if
> I do :
>
>
> ##########################################################################################
> In [31]: master_clump.data.fields Out[31]: ['Density', 'tempContours',
> 'CellMassMsun']
> In [32]: a = master_clump.data['CellMassMsun']
> In [33]: sum(a) Out[33]: 0.98881452898337696
> ##########################################################################################
>
> this will be the total bound mass. Does it automatically take into account
> the thermal energy to determine whether a clump is bound or not ?
>
>
> I would really appreciate if someone could give me a hint.
> Thanks a lot,
>
> Jean-Claude
>
> _______________________________________________
> 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