[yt-users] cloud in cell mesh construction for particle data
Nathan Goldbaum
nathan12343 at gmail.com
Thu Jun 5 18:44:43 PDT 2014
That's right, you can set that via the bbox keyword parameter for
load_particles. I'd urge you to take a look at the docstrings and source
code for load_particles.
On Thu, Jun 5, 2014 at 6:34 PM, Brendan Griffen <brendan.f.griffen at gmail.com
> wrote:
> Thanks very much Nathan. I tried to load in my own data but I think there
> are too many particles or I have to specifically set the domain size.
>
> In this area:
>
> data = {'particle_position_x': pos[:,0],
> 'particle_position_y': pos[:,1],
> 'particle_position_z': pos[:,2],
> 'particle_mass': np.array([mpart]*npart)}
>
> ds = yt.load_particles(data, length_unit=3.08e24, mass_unit=1.9891e36)
> ad = ds.all_data()
> print ad['deposit', 'all_cic']
>
> In [3]: run ytcic.py
> yt : [INFO ] 2014-06-05 21:29:06,183 Parameters: current_time
> = 0.0
> yt : [INFO ] 2014-06-05 21:29:06,183 Parameters: domain_dimensions
> = [2 2 2]
> yt : [INFO ] 2014-06-05 21:29:06,184 Parameters: domain_left_edge
> = [ 0. 0. 0.]
> yt : [INFO ] 2014-06-05 21:29:06,185 Parameters: domain_right_edge
> = [ 1. 1. 1.]
> yt : [INFO ] 2014-06-05 21:29:06,185 Parameters:
> cosmological_simulation = 0.0
> yt : [INFO ] 2014-06-05 21:29:06,188 Allocating for 1.342e+08 particles
> ---------------------------------------------------------------------------
> YTDomainOverflow Traceback (most recent call last)
> /nfs/blank/h4231/bgriffen/data/lib/yt-x86_64/lib/python2.7/site-packages/IPython/utils/py3compat.pyc
> in execfile(fname, *where)
> 202 else:
> 203 filename = fname
> --> 204 __builtin__.execfile(filename, *where)
>
> /nfs/blank/h4231/bgriffen/work/projects/caterpillar/c2ray/cic/ytcic.py in
> <module>()
> 52
> 53 ad = ds.all_data()
> ---> 54 print ad['deposit', 'all_cic']
> 55 slc = yt.SlicePlot(ds, 2, ('deposit', 'all_cic'))
> 56 slc.set_figure_size(4)
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/data_objects/data_containers.pyc
> in __getitem__(self, key)
> 205 Returns a single field. Will add if necessary.
> 206 """
> --> 207 f = self._determine_fields([key])[0]
> 208 if f not in self.field_data and key not in self.field_data:
> 209 if f in self._container_fields:
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/data_objects/data_containers.pyc
> in _determine_fields(self, fields)
> 453 raise YTFieldNotParseable(field)
> 454 ftype, fname = field
> --> 455 finfo = self.pf._get_field_info(ftype, fname)
> 456 else:
> 457 fname = field
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/data_objects/static_output.pyc
> in _get_field_info(self, ftype, fname)
> 445 _last_finfo = None
> 446 def _get_field_info(self, ftype, fname = None):
> --> 447 self.index
> 448 if fname is None:
> 449 ftype, fname = "unknown", ftype
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/data_objects/static_output.pyc
> in index(self)
> 277 raise RuntimeError("You should not instantiate
> Dataset.")
> 278 self._instantiated_index = self._index_class(
> --> 279 self, dataset_type=self.dataset_type)
> 280 # Now we do things that we need an instantiated index
> for
> 281 # ...first off, we create our field_info now.
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/frontends/stream/data_structures.pyc
> in __init__(self, pf, dataset_type)
> 942 def __init__(self, pf, dataset_type = None):
> 943 self.stream_handler = pf.stream_handler
> --> 944 super(StreamParticleIndex, self).__init__(pf, dataset_type)
> 945
> 946 def _setup_data_io(self):
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/particle_geometry_handler.pyc
> in __init__(self, pf, dataset_type)
> 48 self.directory = os.path.dirname(self.index_filename)
> 49 self.float_type = np.float64
> ---> 50 super(ParticleIndex, self).__init__(pf, dataset_type)
> 51
> 52 def _setup_geometry(self):
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/geometry_handler.pyc
> in __init__(self, pf, dataset_type)
> 54
> 55 mylog.debug("Setting up domain geometry.")
> ---> 56 self._setup_geometry()
> 57
> 58 mylog.debug("Initializing data grid data IO")
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/particle_geometry_handler.pyc
> in _setup_geometry(self)
> 52 def _setup_geometry(self):
> 53 mylog.debug("Initializing Particle Geometry Handler.")
> ---> 54 self._initialize_particle_handler()
> 55
> 56
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/particle_geometry_handler.pyc
> in _initialize_particle_handler(self)
> 87 pf.domain_left_edge, pf.domain_right_edge,
> 88 [N, N, N], len(self.data_files))
> ---> 89 self._initialize_indices()
> 90 self.oct_handler.finalize()
> 91 self.max_level = self.oct_handler.max_level
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/particle_geometry_handler.pyc
> in _initialize_indices(self)
> 109 npart = sum(data_file.total_particles.values())
> 110 morton[ind:ind + npart] = \
> --> 111 self.io._initialize_index(data_file, self.regions)
> 112 ind += npart
> 113 morton.sort()
>
> /bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/frontends/stream/io.pyc
> in _initialize_index(self, data_file, regions)
> 144 raise YTDomainOverflow(pos.min(axis=0),
> pos.max(axis=0),
> 145
> data_file.pf.domain_left_edge,
> --> 146
> data_file.pf.domain_right_edge)
> 147 regions.add_data_file(pos, data_file.file_id)
> 148 morton.append(compute_morton(
>
> YTDomainOverflow: Particle bounds [ 0. 0. 0.] and [ 99.99999237
> 99.99999237 99.99999237] exceed domain bounds [ 0. 0. 0.] code_length
> and [ 1. 1. 1.] code_length
>
>
>
> On Thu, Jun 5, 2014 at 8:22 PM, Nathan Goldbaum <nathan12343 at gmail.com>
> wrote:
>
>> Here's a worked out example that does what you're looking for using a
>> fake 1 million particle dataset:
>>
>> http://nbviewer.ipython.org/gist/ngoldbaum/546d37869aafe71cfe38
>>
>> In this notebook I make use of two key yt features: `load_particles`, and
>> `covering_grid`.
>>
>> load_particles creates a "stream" dataset based on in-memory data fed in
>> as a numpy array. This dataset acts just like an on-disk simulation
>> dataset, but doesn't come with the baggage of needing to write a custom
>> frontend to read a specific data format off disk.
>>
>> covering_grid is a way to generate uniform resolution data from an AMR
>> dataset. It acts like a python dictionary where the keys are field names
>> and returns 3D numpy arrays of whatever uniform resolution you specify when
>> you create the covering_grid.
>>
>> Note that if you're using load_particles all of your data needs to live
>> in memory. If your data is too big for that you'll need to write a
>> frontend for your data format or use a memmap to an on-disk file somehow.
>> I'm not an expert on that but others on the list should be able to help
>> out.
>>
>> Hope that gets you well on your way :)
>>
>> -Nathan
>>
>>
>> On Thu, Jun 5, 2014 at 5:04 PM, Desika Narayanan <dnarayan at haverford.edu>
>> wrote:
>>
>>> Hey Brendan,
>>>
>>> A couple of extra tools you might find helpful in conjunction with
>>> Nathan's example of depositing the particles onto an octree are at:
>>>
>>> http://paste.yt-project.org/show/4737/
>>>
>>> Where I load a gadget snapshot, and then recover the coordinates and
>>> width of each cell.
>>>
>>> In response to your last question - the particles are deposited into an
>>> octree grid (so, you'll see that the cell sizes aren't all the same size).
>>> I don't know if depositing onto a regular NxNxN mesh is possible, though
>>> would be interested to hear if so.
>>>
>>> -d
>>>
>>>
>>> On Thu, Jun 5, 2014 at 7:58 PM, Brendan Griffen <
>>> brendan.f.griffen at gmail.com> wrote:
>>>
>>>> Thanks. I'll get the "bleeding edge" version first then try your
>>>> suggestions. Though I want to return the NxNxN array and be able to write
>>>> this mesh to a file. It is *only* using the cic part of yt and it should
>>>> return the mesh to be written? Just wanted to clarify?
>>>>
>>>> Thanks.
>>>> Brendan
>>>>
>>>>
>>>> On Thu, Jun 5, 2014 at 6:49 PM, Nathan Goldbaum <nathan12343 at gmail.com>
>>>> wrote:
>>>>
>>>>>
>>>>>
>>>>>
>>>>> On Thu, Jun 5, 2014 at 3:36 PM, John ZuHone <jzuhone at gmail.com> wrote:
>>>>>
>>>>>> Hi Brendan,
>>>>>>
>>>>>> Which version of yt are you using?
>>>>>>
>>>>>> If you're using 3.0, this is actually fairly easy. If you look in
>>>>>> yt.fields.particle_fields.py, around line 85, you can see how this
>>>>>> is done for the "particle_density" and "particle_mass" fields. Basically
>>>>>> you can call a "deposit" method which takes the particle field quantity you
>>>>>> want deposited and deposits it into cells. The underlying calculation is
>>>>>> done using Cython, so it's fast.
>>>>>>
>>>>>
>>>>> And you shouldn't ever actually need to call these "deposit"
>>>>> functions, since "deposit" is exposed as a field type for all datasets that
>>>>> contain particles.
>>>>>
>>>>> Here is a notebook that does this for Enzo AMR data:
>>>>>
>>>>> http://nbviewer.ipython.org/gist/ngoldbaum/5e19e4e6cc2bf330149c
>>>>>
>>>>> This dataset contains about a million particles and generates a CIC
>>>>> deposition for the whole domain in about 6 seconds from a cold start.
>>>>>
>>>>>
>>>>>>
>>>>>> If you're using 2.x, then you can do the same thing, but it's not as
>>>>>> straightforward. You can see how this works in
>>>>>> yt.data_objects.universal_fields.py, around line 986, where the
>>>>>> "particle_density" field is defined. Basically, it calls CICDeposit_3,
>>>>>> which is in yt.utilities.lib.CICDeposit.pyx.
>>>>>>
>>>>>> Let me know if you need any more clarification.
>>>>>>
>>>>>> Best,
>>>>>>
>>>>>> John Z
>>>>>>
>>>>>> On Jun 5, 2014, at 6:07 PM, Brendan Griffen <
>>>>>> brendan.f.griffen at gmail.com> wrote:
>>>>>>
>>>>>> > Hi,
>>>>>> >
>>>>>> > I was wondering if there were any Cython routines within yt which
>>>>>> takes particle data and converts it into a cloud-in-cell based mesh which
>>>>>> can be written to a file of my choosing.
>>>>>>
>>>>>
>>>>> What sort of mesh were you looking for? yt will internally construct
>>>>> an octree if it is fed particle data. I'm not sure whether this octree can
>>>>> be saved to disk for later analysis.
>>>>>
>>>>> It's also possible to create a uniform resolution covering grid
>>>>> containing field data for a deposited quantity, which can be quite easily
>>>>> saved to disk in a number of ways.
>>>>>
>>>>>
>>>>>> I heard a while ago there was some such functionality but it could be
>>>>>> too far down the yt rabbit hole to be used as a standalone? Is this true? I
>>>>>> have my own Python code for doing it but it just isn't fast enough and
>>>>>> thought I'd ask the yt community if there were any wrapper tools available
>>>>>> to boost the speed.
>>>>>> >
>>>>>> > Thanks.
>>>>>> > Brendan
>>>>>> > _______________________________________________
>>>>>> > yt-users mailing list
>>>>>> > yt-users at lists.spacepope.org
>>>>>> > http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>
>>>>>> _______________________________________________
>>>>>> yt-users mailing list
>>>>>> yt-users at lists.spacepope.org
>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>
>>>>>
>>>>>
>>>>> _______________________________________________
>>>>> yt-users mailing list
>>>>> yt-users at lists.spacepope.org
>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>
>>>>>
>>>>
>>>> _______________________________________________
>>>> yt-users mailing list
>>>> yt-users at lists.spacepope.org
>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>
>>>>
>>>
>>> _______________________________________________
>>> yt-users mailing list
>>> yt-users at lists.spacepope.org
>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>
>>>
>>
>> _______________________________________________
>> yt-users mailing list
>> yt-users at lists.spacepope.org
>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>
>>
>
> _______________________________________________
> yt-users mailing list
> yt-users at lists.spacepope.org
> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.spacepope.org/pipermail/yt-users-spacepope.org/attachments/20140605/e17a7d9b/attachment.html>
More information about the yt-users
mailing list