[Yt-svn] yt: Adding in particle mass function, fixing some issues with le...

[hg at spacepope.org]

hg Repository: yt
details:   yt/rev/0eae5b1c68d5
changeset: 3381:0eae5b1c68d5
user:      Matthew Turk <matthewturk at gmail.com>
date:
Fri Sep 03 16:33:46 2010 -0700
description:
Adding in particle mass function, fixing some issues with left/right edges and
IO for Gadget

diffstat:

 yt/frontends/gadget/data_structures.py |  36 ++++++++++++++++++++++--------------
 yt/frontends/gadget/fields.py          |  17 ++++++++++++++---
 yt/frontends/gadget/io.py              |   2 +-
 3 files changed, 37 insertions(+), 18 deletions(-)

diffs (124 lines):

diff -r a16fa728b1fc -r 0eae5b1c68d5 yt/frontends/gadget/data_structures.py
--- a/yt/frontends/gadget/data_structures.py	Fri Sep 03 15:30:24 2010 -0700
+++ b/yt/frontends/gadget/data_structures.py	Fri Sep 03 16:33:46 2010 -0700
@@ -41,13 +41,13 @@
 
 class GadgetGrid(AMRGridPatch):
     _id_offset = 0
-    def __init__(self, hierarchy, id,dimensions,left_index,
+    def __init__(self, hierarchy, id, dimensions, left_index,
                  level, parent_id, particle_count):
         AMRGridPatch.__init__(self, id, filename = hierarchy.filename,
                               hierarchy = hierarchy)
         self.id = id
         self.ActiveDimensions = dimensions
-        self.start_index = left_index
+        self.start_index = left_index.astype("int64")
         self.Level = level
         self._parent_id = parent_id
         self.Parent = None # Only one parent per grid
@@ -101,24 +101,29 @@
     def _parse_hierarchy(self):
         f = self._handle # shortcut
         npa = na.array
+        DLE = self.parameter_file.domain_left_edge
+        DRE = self.parameter_file.domain_right_edge
+        DW = (DRE - DLE)
         
         self.grid_levels.flat[:] = f['/grid_level'][:].astype("int32")
-        self.grid_left_edge[:]  = f['/grid_left_index'][:]
+        LI = f['/grid_left_index'][:]
+        print LI
         self.grid_dimensions[:] = f['/grid_dimensions'][:]
-        self.grid_right_edge[:] = self.grid_left_edge + self.grid_dimensions 
+        self.grid_left_edge[:]  = (LI * DW + DLE)
+        dxs = 1.0/(2**(self.grid_levels+1)) * DW
+        self.grid_right_edge[:] = self.grid_left_edge \
+                                + dxs *(1 + self.grid_dimensions)
         self.grid_particle_count.flat[:] = f['/grid_particle_count'][:].astype("int32")
         grid_parent_id = f['/grid_parent_id'][:]
         self.max_level = na.max(self.grid_levels)
         
-        args = izip(xrange(self.num_grids), self.grid_levels,
-                    grid_parent_id, self.grid_left_edge,
-                    self.grid_dimensions, self.grid_particle_count)
+        args = izip(xrange(self.num_grids), self.grid_levels.flat,
+                    grid_parent_id, LI,
+                    self.grid_dimensions, self.grid_particle_count.flat)
         self.grids = na.array([self.grid(self,j,d,le,lvl,p,n)
                                for j,lvl,p, le, d, n in args],
                            dtype='object')
         
-        
-        
     def _populate_grid_objects(self):    
         for g in self.grids:
             if g._parent_id >= 0:
@@ -183,9 +188,12 @@
         format = 'Gadget Infrastructure'
         add1 = 'griddded_data_format'
         add2 = 'data_software'
-        fileh = h5py.File(args[0],'r')
-        if add1 in fileh['/'].items():
-            if add2 in fileh['/'+add1].attrs.keys():
-                if fileh['/'+add1].attrs[add2] == format:
-                    return True
+        try:
+            fileh = h5py.File(args[0],'r')
+            if add1 in fileh['/'].items():
+                if add2 in fileh['/'+add1].attrs.keys():
+                    if fileh['/'+add1].attrs[add2] == format:
+                        return True
+        except h5py.h5e.LowLevelIOError:
+            pass
         return False
diff -r a16fa728b1fc -r 0eae5b1c68d5 yt/frontends/gadget/fields.py
--- a/yt/frontends/gadget/fields.py	Fri Sep 03 15:30:24 2010 -0700
+++ b/yt/frontends/gadget/fields.py	Fri Sep 03 16:33:46 2010 -0700
@@ -23,6 +23,9 @@
   along with this program.  If not, see <http://www.gnu.org/licenses/>.
 """
 
+import numpy as na
+
+from yt.funcs import *
 from yt.data_objects.field_info_container import \
     CodeFieldInfoContainer, \
     ValidateParameter, \
@@ -90,13 +93,21 @@
           validators = [ValidateDataField("VEL")],
           units=r"")
 
-add_field("ID", function=lambda a,b: None, take_log=True,
+add_field("id", function=lambda a,b: None, take_log=True,
           validators = [ValidateDataField("ID")],
           units=r"")
 
-add_field("MASS", function=lambda a,b: None, take_log=True,
-          validators = [ValidateDataField("MASS")],
+add_field("mass", function=lambda a,b: None, take_log=True,
+          validators = [ValidateDataField("mass")],
           units=r"\rm{g}")
+def _particle_mass(field, data):
+    return data["mass"]/just_one(data["CellVolume"])
+def _convert_particle_mass(data):
+    return 1.0
+add_field("particle_mass", function=_particle_mass, take_log=True,
+          convert_function=_convert_particle_mass,
+          validators = [ValidateSpatial(0)],
+          units=r"\mathrm{g}\/\mathrm{cm}^{-3}")
 
 add_field("U", function=lambda a,b: None, take_log=True,
           validators = [ValidateDataField("U")],
diff -r a16fa728b1fc -r 0eae5b1c68d5 yt/frontends/gadget/io.py
--- a/yt/frontends/gadget/io.py	Fri Sep 03 15:30:24 2010 -0700
+++ b/yt/frontends/gadget/io.py	Fri Sep 03 16:33:46 2010 -0700
@@ -38,7 +38,7 @@
             address = '/data/grid_%010i/particles/%s/%s' % (grid.id, ptype, field)
             data.append(fh[address][:])
         if len(data) > 0:
-            data = na.vstack(data)
+            data = na.concatenate(data)
         fh.close()
         return na.array(data)
     def _read_field_names(self,grid):