[Yt-svn] yt-commit r593 - branches/yt-generalization/yt/lagos

joishi at wrangler.dreamhost.com joishi at wrangler.dreamhost.com
Sat Jun 21 13:46:14 PDT 2008


Author: joishi
Date: Sat Jun 21 13:46:13 2008
New Revision: 593
URL: http://yt.spacepope.org/changeset/593

Log:
* a first cut at Orion reading. Doesn't work, but doesn't interrupt Enzo tests.


Modified:
   branches/yt-generalization/yt/lagos/DataReadingFuncs.py
   branches/yt-generalization/yt/lagos/HierarchyType.py
   branches/yt-generalization/yt/lagos/OutputTypes.py

Modified: branches/yt-generalization/yt/lagos/DataReadingFuncs.py
==============================================================================
--- branches/yt-generalization/yt/lagos/DataReadingFuncs.py	(original)
+++ branches/yt-generalization/yt/lagos/DataReadingFuncs.py	Sat Jun 21 13:46:13 2008
@@ -142,3 +142,14 @@
 
 def getExceptionHDF5():
     return (exceptions.KeyError, HDF5LightReader.ReadingError)
+
+def readDataNative():
+    pass
+
+def readAllDataNative():
+    pass
+
+def readDataSliceNative():
+    """wishful thinking?
+    """
+    pass

Modified: branches/yt-generalization/yt/lagos/HierarchyType.py
==============================================================================
--- branches/yt-generalization/yt/lagos/HierarchyType.py	(original)
+++ branches/yt-generalization/yt/lagos/HierarchyType.py	Sat Jun 21 13:46:13 2008
@@ -32,7 +32,8 @@
 _data_style_funcs = \
    { 4: (readDataHDF4, readAllDataHDF4, getFieldsHDF4, readDataSliceHDF4, getExceptionHDF4), \
      5: (readDataHDF5, readAllDataHDF5, getFieldsHDF5, readDataSliceHDF5, getExceptionHDF5), \
-     6: (readDataPacked, readAllDataPacked, getFieldsPacked, readDataSlicePacked, getExceptionHDF5) \
+     6: (readDataPacked, readAllDataPacked, getFieldsPacked, readDataSlicePacked, getExceptionHDF5), \
+     7: (readDataNative, readAllDataNative, readDataSliceNative) \
    }
 
 class AMRHierarchy:
@@ -796,3 +797,94 @@
         rs += "(%s)\s*" % (scanf_regex[t])
     rs +="$"
     return re.compile(rs,re.M)
+
+class OrionHierarchy(AMRHierarchy):
+    def __init__(self,filename):
+        self.readGlobalHeader(filename)
+        
+    def readGlobalHeader(self,filename):
+        """
+        read the global header file for an Orion plotfile output.
+        
+        
+        """
+
+        # domain_pattern = r"\(\((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\)\) \(\((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\)\)"
+        # domain_re = re.compile(domain_pattern)
+
+        counter = 0
+        # CHANGE FOLLOWING: __header_file SHOULD NOT BE A CLASS MEMBER
+        self.__header_file = open(filename,'r')
+        self.__global_header_lines = self.__header_file.readlines()
+
+        # parse the file
+        self.orion_version = self.__global_header_lines[0].rstrip()
+        self.n_fields      = int(self.__global_header_lines[1])
+
+        self.fields = []
+        counter = self.n_fields+2
+        for line in self.__global_header_lines[2:counter]:
+            self.fields.append(line.rstrip())
+
+        self.dimension = int(self.__global_header_lines[counter])
+        if self.dimension != 3:
+            raise RunTimeError("Orion must be in 3D to use yt.")
+        counter += 1
+        self.Time = float(self.__global_header_lines[counter])
+        counter += 1
+        self.finest_grid_level = int(self.__global_header_lines[counter])
+        self.n_levels = self.finest_grid_level + 1
+        counter += 1
+        self.domainLeftEdge_unnecessary = na.array(map(float,self.__global_header_lines[counter].split()))
+        counter += 1
+        self.domainRightEdge_unnecessary = na.array(map(float,self.__global_header_lines[counter].split()))
+        counter += 1
+        self.refinementFactor_unnecessary = na.array(map(int,self.__global_header_lines[counter].split()))
+        counter += 1
+        self.globalIndexSpace_unnecessary = self.__global_header_lines[counter]
+        #domain_re.search(self.__global_header_lines[counter]).groups()
+        counter += 1
+        self.timestepsPerLevel_unnecessary = self.__global_header_lines[counter]
+        counter += 1
+        self.dx = na.zeros((self.n_levels,3))
+        for i,line in enumerate(self.__global_header_lines[counter:counter+self.n_levels]):
+            self.dx[i] = na.array(map(int,line.split()))
+        counter += self.n_levels
+        self.geometry = int(self.__global_header_lines[counter])
+        if self.geometry != 0:
+            raise RunTimeError("yt only supports cartesian coordinates.")
+        counter += 1
+
+        # this is just to debug. eventually it should go away.
+        linebreak = int(self.__global_header_lines[counter])
+        if linebreak != 0:
+            raise RunTimeError("INTERNAL ERROR! This should be a zero. the reader function be broke, J.S.")
+        counter += 1
+
+        # each level is one group with ngrids on it. each grid has 3 lines of 2 reals
+        self.levels = []
+        grid_counter = 0
+        for level in range(0,self.n_levels):
+            tmp = self.__global_header_lines[counter].split()
+            # should this be grid_time or level_time??
+            lev,ngrids,grid_time = int(tmp[0]),int(tmp[1]),float(tmp[2])
+            counter += 1
+            nsteps = int(self.__global_header_lines[counter])
+            counter += 1
+            self.levels.append(OrionLevel(lev,ngrids))
+            for grid in range(0,ngrids):
+                xlo,xhi = map(float,self.__global_header_lines[counter].split())
+                counter+=1
+                ylo,yhi = map(float,self.__global_header_lines[counter].split())
+                counter+=1
+                zlo,zhi = map(float,self.__global_header_lines[counter].split())
+                counter+=1
+                lo = na.array([xlo,ylo,zlo])
+                hi = na.array([xhi,yhi,zhi])
+                self.levels[-1].grids.append(OrionGrid(lo,hi,grid_counter))
+                grid_counter += 1 # this is global, and shouldn't be reset
+                                  # for each level
+            self.levels[-1]._fileprefix = self.__global_header_lines[counter]
+            counter+=1
+
+        self.__header_file.close()

Modified: branches/yt-generalization/yt/lagos/OutputTypes.py
==============================================================================
--- branches/yt-generalization/yt/lagos/OutputTypes.py	(original)
+++ branches/yt-generalization/yt/lagos/OutputTypes.py	Sat Jun 21 13:46:13 2008
@@ -6,7 +6,7 @@
 @organization: U{KIPAC<http://www-group.slac.stanford.edu/KIPAC/>}
 @contact: U{mturk at slac.stanford.edu<mailto:mturk at slac.stanford.edu>}
 @license:
-  Copyright (C) 2007 Matthew Turk.  All Rights Reserved.
+  Copyright (C) 2007-2008 Matthew Turk, J. S. Oishi.  All Rights Reserved.
 
   This file is part of yt.
 
@@ -269,3 +269,124 @@
         k["aye"]  = (1.0 + self.parameters["CosmologyInitialRedshift"]) / \
                (1.0 + self.parameters["CosmologyCurrentRedshift"])
         return k
+
+class OrionStaticOutput(StaticOutput):
+    """
+    This class is a stripped down class that simply reads and parses, without
+    looking at the Orion hierarchy.
+
+    @todo: 
+
+    @param filename: The filename of the parameterfile we want to load
+    @type filename: String
+    """
+    _hierarchy_class = OrionHierarchy
+
+    def __init__(self, plotname, paramFilename='inputs',data_style='Native'):
+        """need to override for Orion file structure.
+
+        the paramfile is usually called "inputs"
+        plotname here will be a directory name
+        as per BoxLib, data_style will be one of
+          Native
+          IEEE
+          ASCII
+
+        """
+        self.data_style = data_style
+        self.basename = os.path.basename(plotname)
+        # this will be the directory ENCLOSING the pltNNNN directory
+        self.directory = os.path.dirname(plotname)
+        self.parameter_filename = self.directory+paramFilename
+        self.fullpath = os.path.abspath(self.directory)
+        if len(self.directory) == 0:
+            self.directory = "."
+        self.conversion_factors = {}
+        self.parameters = {}
+        self._parse_parameter_file()
+        self._set_units()
+        
+        # These should maybe not be hardcoded?
+        self.parameters["HydroMethod"] = 0 # always PPM DE
+        self.parameters["InitialTime"] = 0. # FIX ME!!!
+        
+    def _parse_parameter_file(self):
+        """
+        Parses the parameter file and establishes the various
+        dictionaries.
+        """
+        # Let's read the file
+        self.parameters["CurrentTimeIdentifier"] = \
+            int(os.stat(self.parameter_filename)[ST_CTIME])
+        lines = open(self.parameter_filename).readlines()
+        for lineI, line in enumerate(lines):
+            if line.find("#") >= 1: # Keep the commented lines...
+                line=line[:line.find("#")]
+            line=line.strip().rstrip()
+            if len(line) < 2 or line.find("#") == 0: # ...but skip comments
+                continue
+            try:
+                param, vals = map(strip,map(rstrip,line.split("=")))
+            except ValueError:
+                mylog.error("ValueError: '%s'", line)
+            if orion2enzoDict.has_key(param):
+                paramName = orion2enzoDict[param]
+                t = map(parameterDict[paramName], vals.split())
+                if len(t) == 1:
+                    self.parameters[paramName] = t[0]
+                else:
+                    self.parameters[paramName] = t
+                if paramName.endswith("Units") and not paramName.startswith("Temperature"):
+                    dataType = paramName[:-5]
+                    self.conversion_factors[dataType] = self.parameters[paramName]
+            # elif param.startswith("#DataCGS"):
+#                 # Assume of the form: #DataCGSConversionFactor[7] = 2.38599e-26 g/cm^3
+#                 if lines[lineI-1].find("Label") >= 0:
+#                     kk = lineI-1
+#                 elif lines[lineI-2].find("Label") >= 0:
+#                     kk = lineI-2
+#                 dataType = lines[kk].split("=")[-1].rstrip().strip()
+#                 convFactor = float(line.split("=")[-1].split()[0])
+#                 self.conversion_factors[dataType] = convFactor
+#             elif param.startswith("#CGSConversionFactor"):
+#                 dataType = param[20:].rstrip()
+#                 convFactor = float(line.split("=")[-1])
+#                 self.conversion_factors[dataType] = convFactor
+            elif param.startswith("geometry.prob_hi"):
+                self.parameters["DomainRightEdge"] = \
+                    na.array([float(i) for i in vals.split()])
+            elif param.startswith("geometry.prob_lo"):
+                self.parameters["DomainLeftEdge"] = \
+                    na.array([float(i) for i in vals.split()])
+                
+    def _set_units(self):
+        """
+        Generates the conversion to various physical _units based on the parameter file
+        """
+        self.units = {}
+        self.time_units = {}
+        if len(self.parameters) == 0:
+            self._parse_parameter_file()
+        if self["ComovingCoordinates"]:
+            raise RuntimeError("Orion doesn't do cosmology.")
+        elif self.has_key("LengthUnits"):
+            raise RuntimeError("Units not yet implemented.")
+            #self._setup_getunits_units()
+        else:
+            self._setup_nounits_units()
+        self.time_units['1'] = 1
+        self.units['1'] = 1
+        seconds = self["Time"]
+        self.time_units['years'] = seconds / (365*3600*24.0)
+        self.time_units['days']  = seconds / (3600*24.0)
+
+    def _setup_nounits_units(self):
+        z = 0
+        box_proper = ytcfg.getfloat("lagos","nounitslength")
+        self.units['aye'] = 1.0
+        mylog.warning("No length units.  Setting 1.0 = %0.3e proper Mpc.", box_proper)
+        if not self.has_key("TimeUnits"):
+            mylog.warning("No time units.  Setting 1.0 = 1 second.")
+            self.conversion_factors["Time"] = 1.0
+        for unit in mpc_conversion.keys():
+            self.units[unit] = mpc_conversion[unit] * box_proper



More information about the yt-svn mailing list