[yt-svn] commit/yt: 3 new changesets
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Wed May 15 12:15:28 PDT 2013
3 new commits in yt:
https://bitbucket.org/yt_analysis/yt/commits/b2cb6f33ba48/
Changeset: b2cb6f33ba48
Branch: yt
User: MatthewTurk
Date: 2013-05-15 19:05:53
Summary: Copying from Orion into new Boxlib frontend. Will merge Castro, Maestro and Nyx.
Affected #: 8 files
diff -r a03fc16b5457aa3deba49a11277fbb41ffebd09e -r b2cb6f33ba486efad0f9a62f2c523900c899a65e yt/frontends/boxlib/__init__.py
--- /dev/null
+++ b/yt/frontends/boxlib/__init__.py
@@ -0,0 +1,29 @@
+"""
+API for yt.frontends.orion
+
+Author: Matthew Turk <matthewturk at gmail.com>
+Affiliation: UCSD
+Author: J.S. Oishi <jsoishi at gmail.com>
+Affiliation: KIPAC/SLAC/Stanford
+Author: Britton Smith <brittonsmith at gmail.com>
+Affiliation: MSU
+Homepage: http://yt-project.org/
+License:
+ Copyright (C) 2010-2011 Matthew Turk. All Rights Reserved.
+
+ This file is part of yt.
+
+ yt is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+"""
diff -r a03fc16b5457aa3deba49a11277fbb41ffebd09e -r b2cb6f33ba486efad0f9a62f2c523900c899a65e yt/frontends/boxlib/api.py
--- /dev/null
+++ b/yt/frontends/boxlib/api.py
@@ -0,0 +1,41 @@
+"""
+API for yt.frontends.orion
+
+Author: Matthew Turk <matthewturk at gmail.com>
+Affiliation: UCSD
+Author: J.S. Oishi <jsoishi at gmail.com>
+Affiliation: KIPAC/SLAC/Stanford
+Author: Britton Smith <brittonsmith at gmail.com>
+Affiliation: MSU
+Homepage: http://yt-project.org/
+License:
+ Copyright (C) 2010-2011 Matthew Turk. All Rights Reserved.
+
+ This file is part of yt.
+
+ yt is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+"""
+
+from .data_structures import \
+ OrionGrid, \
+ OrionHierarchy, \
+ OrionStaticOutput
+
+from .fields import \
+ OrionFieldInfo, \
+ add_orion_field
+
+from .io import \
+ IOHandlerNative
diff -r a03fc16b5457aa3deba49a11277fbb41ffebd09e -r b2cb6f33ba486efad0f9a62f2c523900c899a65e yt/frontends/boxlib/data_structures.py
--- /dev/null
+++ b/yt/frontends/boxlib/data_structures.py
@@ -0,0 +1,645 @@
+"""
+Data structures for Orion.
+
+Author: J. S. Oishi <jsoishi at gmail.com>
+Affiliation: KIPAC/SLAC/Stanford
+Homepage: http://yt-project.org/
+License:
+ Copyright (C) 2008-2011 J. S. Oishi. All Rights Reserved.
+
+ This file is part of yt.
+
+ yt is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+"""
+
+import os
+import re
+import weakref
+
+from collections import defaultdict
+from string import strip, rstrip
+from stat import ST_CTIME
+
+import numpy as np
+
+from yt.funcs import *
+from yt.data_objects.field_info_container import FieldInfoContainer, NullFunc
+from yt.data_objects.grid_patch import AMRGridPatch
+from yt.data_objects.hierarchy import AMRHierarchy
+from yt.data_objects.static_output import StaticOutput
+from yt.utilities.definitions import \
+ mpc_conversion, sec_conversion
+from yt.utilities.parallel_tools.parallel_analysis_interface import \
+ parallel_root_only
+
+from .definitions import \
+ orion2enzoDict, \
+ parameterDict, \
+ yt2orionFieldsDict, \
+ orion_FAB_header_pattern
+from .fields import \
+ OrionFieldInfo, \
+ add_orion_field, \
+ KnownOrionFields
+
+
+class OrionGrid(AMRGridPatch):
+ _id_offset = 0
+ def __init__(self, LeftEdge, RightEdge, index, level, filename, offset,
+ dimensions, start, stop, paranoia=False, **kwargs):
+ AMRGridPatch.__init__(self, index, **kwargs)
+ self.filename = filename
+ self._offset = offset
+ self._paranoid = paranoia
+
+ # should error check this
+ self.ActiveDimensions = (dimensions.copy()).astype('int32')#.transpose()
+ self.start_index = start.copy()#.transpose()
+ self.stop_index = stop.copy()#.transpose()
+ self.LeftEdge = LeftEdge.copy()
+ self.RightEdge = RightEdge.copy()
+ self.index = index
+ self.Level = level
+
+ def get_global_startindex(self):
+ return self.start_index
+
+ def _prepare_grid(self):
+ """
+ Copies all the appropriate attributes from the hierarchy
+ """
+ # This is definitely the slowest part of generating the hierarchy
+ # Now we give it pointers to all of its attributes
+ # Note that to keep in line with Enzo, we have broken PEP-8
+ h = self.hierarchy # cache it
+ #self.StartIndices = h.gridStartIndices[self.id]
+ #self.EndIndices = h.gridEndIndices[self.id]
+ h.grid_levels[self.id,0] = self.Level
+ h.grid_left_edge[self.id,:] = self.LeftEdge[:]
+ h.grid_right_edge[self.id,:] = self.RightEdge[:]
+ #self.Time = h.gridTimes[self.id,0]
+ self.NumberOfParticles = 0 # these will be read in later
+ self.field_indexes = h.field_indexes
+ self.Children = h.gridTree[self.id]
+ pIDs = h.gridReverseTree[self.id]
+ if len(pIDs) > 0:
+ self.Parent = [weakref.proxy(h.grids[pID]) for pID in pIDs]
+ else:
+ self.Parent = None
+
+ def _setup_dx(self):
+ # has already been read in and stored in hierarchy
+ dx = self.hierarchy.grid_dxs[self.index][0]
+ dy = self.hierarchy.grid_dys[self.index][0]
+ dz = self.hierarchy.grid_dzs[self.index][0]
+ self.dds = np.array([dx, dy, dz])
+ self.field_data['dx'], self.field_data['dy'], self.field_data['dz'] = self.dds
+
+ def __repr__(self):
+ return "OrionGrid_%04i" % (self.id)
+
+class OrionHierarchy(AMRHierarchy):
+ grid = OrionGrid
+ def __init__(self, pf, data_style='orion_native'):
+ self.field_indexes = {}
+ self.parameter_file = weakref.proxy(pf)
+ header_filename = os.path.join(pf.fullplotdir,'Header')
+ self.directory = pf.fullpath
+ self.data_style = data_style
+
+ self.readGlobalHeader(header_filename,self.parameter_file.paranoid_read) # also sets up the grid objects
+ self.__cache_endianness(self.levels[-1].grids[-1])
+ AMRHierarchy.__init__(self,pf, self.data_style)
+ self._populate_hierarchy()
+ self._read_particles()
+
+ def _read_particles(self):
+ """
+ reads in particles and assigns them to grids. Will search for
+ Star particles, then sink particles if no star particle file
+ is found, and finally will simply note that no particles are
+ found if neither works. To add a new Orion particle type,
+ simply add it to the if/elif/else block.
+
+ """
+ self.grid_particle_count = np.zeros(len(self.grids))
+
+ for particle_filename in ["StarParticles", "SinkParticles"]:
+ fn = os.path.join(self.pf.fullplotdir, particle_filename)
+ if os.path.exists(fn): self._read_particle_file(fn)
+
+ def _read_particle_file(self, fn):
+ """actually reads the orion particle data file itself.
+
+ """
+ if not os.path.exists(fn): return
+ with open(fn, 'r') as f:
+ lines = f.readlines()
+ self.num_stars = int(lines[0].strip()[0])
+ for line in lines[1:]:
+ particle_position_x = float(line.split(' ')[1])
+ particle_position_y = float(line.split(' ')[2])
+ particle_position_z = float(line.split(' ')[3])
+ coord = [particle_position_x, particle_position_y, particle_position_z]
+ # for each particle, determine which grids contain it
+ # copied from object_finding_mixin.py
+ mask=np.ones(self.num_grids)
+ for i in xrange(len(coord)):
+ np.choose(np.greater(self.grid_left_edge[:,i],coord[i]), (mask,0), mask)
+ np.choose(np.greater(self.grid_right_edge[:,i],coord[i]), (0,mask), mask)
+ ind = np.where(mask == 1)
+ selected_grids = self.grids[ind]
+ # in orion, particles always live on the finest level.
+ # so, we want to assign the particle to the finest of
+ # the grids we just found
+ if len(selected_grids) != 0:
+ grid = sorted(selected_grids, key=lambda grid: grid.Level)[-1]
+ ind = np.where(self.grids == grid)[0][0]
+ self.grid_particle_count[ind] += 1
+ self.grids[ind].NumberOfParticles += 1
+ return True
+
+ def readGlobalHeader(self,filename,paranoid_read):
+ """
+ read the global header file for an Orion plotfile output.
+ """
+ counter = 0
+ header_file = open(filename,'r')
+ self.__global_header_lines = header_file.readlines()
+
+ # parse the file
+ self.orion_version = self.__global_header_lines[0].rstrip()
+ self.n_fields = int(self.__global_header_lines[1])
+
+ counter = self.n_fields+2
+ self.field_list = []
+ for i,line in enumerate(self.__global_header_lines[2:counter]):
+ self.field_list.append(line.rstrip())
+
+ # this is unused...eliminate it?
+ #for f in self.field_indexes:
+ # self.field_list.append(orion2ytFieldsDict.get(f,f))
+
+ 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
+ # quantities with _unnecessary are also stored in the inputs
+ # file and are not needed. they are read in and stored in
+ # case in the future we want to enable a "backwards" way of
+ # taking the data out of the Header file and using it to fill
+ # in in the case of a missing inputs file
+ self.domainLeftEdge_unnecessary = np.array(map(float,self.__global_header_lines[counter].split()))
+ counter += 1
+ self.domainRightEdge_unnecessary = np.array(map(float,self.__global_header_lines[counter].split()))
+ counter += 1
+ self.refinementFactor_unnecessary = self.__global_header_lines[counter].split() #np.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 = np.zeros((self.n_levels,3))
+ for i,line in enumerate(self.__global_header_lines[counter:counter+self.n_levels]):
+ self.dx[i] = np.array(map(float,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.")
+ counter += 1
+
+ # each level is one group with ngrids on it. each grid has 3 lines of 2 reals
+ self.levels = []
+ grid_counter = 0
+ file_finder_pattern = r"FabOnDisk: (\w+_D_[0-9]{4}) (\d+)\n"
+ re_file_finder = re.compile(file_finder_pattern)
+ dim_finder_pattern = r"\(\((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\) \(\d+,\d+,\d+\)\)\n"
+ re_dim_finder = re.compile(dim_finder_pattern)
+ data_files_pattern = r"Level_[\d]/"
+ data_files_finder = re.compile(data_files_pattern)
+
+ 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))
+ # open level header, extract file names and offsets for
+ # each grid
+ # read slightly out of order here: at the end of the lo,hi
+ # pairs for x,y,z is a *list* of files types in the Level
+ # directory. each type has Header and a number of data
+ # files (one per processor)
+ tmp_offset = counter + 3*ngrids
+ nfiles = 0
+ key_off = 0
+ files = {} # dict(map(lambda a: (a,[]),self.field_list))
+ offsets = {} # dict(map(lambda a: (a,[]),self.field_list))
+ while nfiles+tmp_offset < len(self.__global_header_lines) and data_files_finder.match(self.__global_header_lines[nfiles+tmp_offset]):
+ filen = os.path.join(self.parameter_file.fullplotdir, \
+ self.__global_header_lines[nfiles+tmp_offset].strip())
+ # open each "_H" header file, and get the number of
+ # components within it
+ level_header_file = open(filen+'_H','r').read()
+ start_stop_index = re_dim_finder.findall(level_header_file) # just take the last one
+ grid_file_offset = re_file_finder.findall(level_header_file)
+ ncomp_this_file = int(level_header_file.split('\n')[2])
+ for i in range(ncomp_this_file):
+ key = self.field_list[i+key_off]
+ f,o = zip(*grid_file_offset)
+ files[key] = f
+ offsets[key] = o
+ self.field_indexes[key] = i
+ key_off += ncomp_this_file
+ nfiles += 1
+ # convert dict of lists to list of dicts
+ fn = []
+ off = []
+ lead_path = os.path.join(self.parameter_file.fullplotdir,'Level_%i'%level)
+ for i in range(ngrids):
+ fi = [os.path.join(lead_path,files[key][i]) for key in self.field_list]
+ of = [int(offsets[key][i]) for key in self.field_list]
+ fn.append(dict(zip(self.field_list,fi)))
+ off.append(dict(zip(self.field_list,of)))
+
+ for grid in range(0,ngrids):
+ gfn = fn[grid] # filename of file containing this grid
+ gfo = off[grid] # offset within that file
+ 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 = np.array([xlo,ylo,zlo])
+ hi = np.array([xhi,yhi,zhi])
+ dims,start,stop = self.__calculate_grid_dimensions(start_stop_index[grid])
+ self.levels[-1].grids.append(self.grid(lo,hi,grid_counter,level,gfn, gfo, dims,start,stop,paranoia=paranoid_read,hierarchy=self))
+ grid_counter += 1 # this is global, and shouldn't be reset
+ # for each level
+
+ # already read the filenames above...
+ counter+=nfiles
+ self.num_grids = grid_counter
+ self.float_type = 'float64'
+
+ self.maxLevel = self.n_levels - 1
+ self.max_level = self.n_levels - 1
+ header_file.close()
+
+ def __cache_endianness(self,test_grid):
+ """
+ Cache the endianness and bytes perreal of the grids by using a
+ test grid and assuming that all grids have the same
+ endianness. This is a pretty safe assumption since Orion uses
+ one file per processor, and if you're running on a cluster
+ with different endian processors, then you're on your own!
+ """
+ # open the test file & grab the header
+ inFile = open(os.path.expanduser(test_grid.filename[self.field_list[0]]),'rb')
+ header = inFile.readline()
+ inFile.close()
+ header.strip()
+
+ # parse it. the patter is in OrionDefs.py
+ headerRe = re.compile(orion_FAB_header_pattern)
+ bytesPerReal,endian,start,stop,centerType,nComponents = headerRe.search(header).groups()
+ self._bytesPerReal = int(bytesPerReal)
+ if self._bytesPerReal == int(endian[0]):
+ dtype = '<'
+ elif self._bytesPerReal == int(endian[-1]):
+ dtype = '>'
+ else:
+ raise ValueError("FAB header is neither big nor little endian. Perhaps the file is corrupt?")
+
+ dtype += ('f%i' % self._bytesPerReal) # always a floating point
+ self._dtype = dtype
+
+ def __calculate_grid_dimensions(self,start_stop):
+ start = np.array(map(int,start_stop[0].split(',')))
+ stop = np.array(map(int,start_stop[1].split(',')))
+ dimension = stop - start + 1
+ return dimension,start,stop
+
+ def _populate_grid_objects(self):
+ mylog.debug("Creating grid objects")
+ self.grids = np.concatenate([level.grids for level in self.levels])
+ self.grid_levels = np.concatenate([level.ngrids*[level.level] for level in self.levels])
+ self.grid_levels = np.array(self.grid_levels.reshape((self.num_grids,1)),dtype='int32')
+ grid_dcs = np.concatenate([level.ngrids*[self.dx[level.level]] for level in self.levels],axis=0)
+ self.grid_dxs = grid_dcs[:,0].reshape((self.num_grids,1))
+ self.grid_dys = grid_dcs[:,1].reshape((self.num_grids,1))
+ self.grid_dzs = grid_dcs[:,2].reshape((self.num_grids,1))
+ left_edges = []
+ right_edges = []
+ dims = []
+ for level in self.levels:
+ left_edges += [g.LeftEdge for g in level.grids]
+ right_edges += [g.RightEdge for g in level.grids]
+ dims += [g.ActiveDimensions for g in level.grids]
+ self.grid_left_edge = np.array(left_edges)
+ self.grid_right_edge = np.array(right_edges)
+ self.grid_dimensions = np.array(dims)
+ self.gridReverseTree = [] * self.num_grids
+ self.gridReverseTree = [ [] for i in range(self.num_grids)]
+ self.gridTree = [ [] for i in range(self.num_grids)]
+ mylog.debug("Done creating grid objects")
+
+ def _populate_hierarchy(self):
+ self.__setup_grid_tree()
+ #self._setup_grid_corners()
+ for i, grid in enumerate(self.grids):
+ if (i%1e4) == 0: mylog.debug("Prepared % 7i / % 7i grids", i, self.num_grids)
+ grid._prepare_grid()
+ grid._setup_dx()
+
+ def __setup_grid_tree(self):
+ for i, grid in enumerate(self.grids):
+ children = self._get_grid_children(grid)
+ for child in children:
+ self.gridReverseTree[child.id].append(i)
+ self.gridTree[i].append(weakref.proxy(child))
+
+ def _setup_classes(self):
+ dd = self._get_data_reader_dict()
+ dd["field_indexes"] = self.field_indexes
+ AMRHierarchy._setup_classes(self, dd)
+ #self._add_object_class('grid', "OrionGrid", OrionGridBase, dd)
+ self.object_types.sort()
+
+ def _get_grid_children(self, grid):
+ mask = np.zeros(self.num_grids, dtype='bool')
+ grids, grid_ind = self.get_box_grids(grid.LeftEdge, grid.RightEdge)
+ mask[grid_ind] = True
+ mask = np.logical_and(mask, (self.grid_levels == (grid.Level+1)).flat)
+ return self.grids[mask]
+
+ def _count_grids(self):
+ """this is already provided in
+
+ """
+ pass
+
+ def _initialize_grid_arrays(self):
+ mylog.debug("Allocating arrays for %s grids", self.num_grids)
+ self.grid_dimensions = np.ones((self.num_grids,3), 'int32')
+ self.grid_left_edge = np.zeros((self.num_grids,3), self.float_type)
+ self.grid_right_edge = np.ones((self.num_grids,3), self.float_type)
+ self.grid_levels = np.zeros((self.num_grids,1), 'int32')
+ self.grid_particle_count = np.zeros((self.num_grids,1), 'int32')
+
+ def _parse_hierarchy(self):
+ pass
+
+ def _detect_fields(self):
+ pass
+
+ def _setup_derived_fields(self):
+ pass
+
+ def _initialize_state_variables(self):
+ """override to not re-initialize num_grids in AMRHierarchy.__init__
+
+ """
+ self._parallel_locking = False
+ self._data_file = None
+ self._data_mode = None
+ self._max_locations = {}
+
+class OrionLevel:
+ def __init__(self,level,ngrids):
+ self.level = level
+ self.ngrids = ngrids
+ self.grids = []
+
+
+class OrionStaticOutput(StaticOutput):
+ """
+ This class is a stripped down class that simply reads and parses
+ *filename*, without looking at the Orion hierarchy.
+ """
+ _hierarchy_class = OrionHierarchy
+ _fieldinfo_fallback = OrionFieldInfo
+ _fieldinfo_known = KnownOrionFields
+
+ def __init__(self, plotname, paramFilename=None, fparamFilename=None,
+ data_style='orion_native', paranoia=False,
+ storage_filename = None):
+ """
+ The paramfile is usually called "inputs"
+ and there may be a fortran inputs file usually called "probin"
+ plotname here will be a directory name
+ as per BoxLib, data_style will be Native (implemented here), IEEE (not
+ yet implemented) or ASCII (not yet implemented.)
+ """
+ self.storage_filename = storage_filename
+ self.paranoid_read = paranoia
+ self.parameter_filename = paramFilename
+ self.fparameter_filename = fparamFilename
+ self.__ipfn = paramFilename
+
+ self.fparameters = {}
+
+ StaticOutput.__init__(self, plotname.rstrip("/"),
+ data_style='orion_native')
+
+ # These should maybe not be hardcoded?
+ self.parameters["HydroMethod"] = 'orion' # always PPM DE
+ self.parameters["Time"] = 1. # default unit is 1...
+ self.parameters["DualEnergyFormalism"] = 0 # always off.
+ self.parameters["EOSType"] = -1 # default
+
+ if self.fparameters.has_key("mu"):
+ self.parameters["mu"] = self.fparameters["mu"]
+
+ def _localize(self, f, default):
+ if f is None:
+ return os.path.join(self.directory, default)
+ return f
+
+ @classmethod
+ def _is_valid(cls, *args, **kwargs):
+ # fill our args
+ pname = args[0].rstrip("/")
+ dn = os.path.dirname(pname)
+ if len(args) > 1: kwargs['paramFilename'] = args[1]
+ pfname = kwargs.get("paramFilename", os.path.join(dn, "inputs"))
+
+ # We check for the job_info file's existence because this is currently
+ # what distinguishes Orion data from MAESTRO data.
+ pfn = os.path.join(pfname)
+ if not os.path.exists(pfn): return False
+ castro = any(("castro." in line for line in open(pfn)))
+ nyx = any(("nyx." in line for line in open(pfn)))
+ maestro = os.path.exists(os.path.join(pname, "job_info"))
+ really_orion = any(("geometry.prob_lo" in line for line in open(pfn)))
+ orion = (not castro) and (not maestro) and (not nyx) and really_orion
+ return orion
+
+ def _parse_parameter_file(self):
+ """
+ Parses the parameter file and establishes the various
+ dictionaries.
+ """
+ self.fullplotdir = os.path.abspath(self.parameter_filename)
+ self._parse_header_file()
+ self.parameter_filename = self._localize(
+ self.__ipfn, 'inputs')
+ self.fparameter_filename = self._localize(
+ self.fparameter_filename, 'probin')
+ if os.path.isfile(self.fparameter_filename):
+ self._parse_fparameter_file()
+ for param in self.fparameters:
+ if orion2enzoDict.has_key(param):
+ self.parameters[orion2enzoDict[param]]=self.fparameters[param]
+ # Let's read the file
+ self.unique_identifier = \
+ 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)
+ continue
+ if orion2enzoDict.has_key(param):
+ paramName = orion2enzoDict[param]
+ t = map(parameterDict[paramName], vals.split())
+ if len(t) == 1:
+ self.parameters[paramName] = t[0]
+ else:
+ if paramName == "RefineBy":
+ self.parameters[paramName] = t[0]
+ else:
+ self.parameters[paramName] = t
+
+ elif param.startswith("geometry.prob_hi"):
+ self.domain_right_edge = \
+ np.array([float(i) for i in vals.split()])
+ elif param.startswith("geometry.prob_lo"):
+ self.domain_left_edge = \
+ np.array([float(i) for i in vals.split()])
+ elif param.startswith("Prob.lo_bc"):
+ self.periodicity = ensure_tuple([i == 0 for i in vals])
+
+ self.parameters["TopGridRank"] = len(self.parameters["TopGridDimensions"])
+ self.dimensionality = self.parameters["TopGridRank"]
+ self.domain_dimensions = np.array(self.parameters["TopGridDimensions"],dtype='int32')
+ self.refine_by = self.parameters["RefineBy"]
+
+ if self.parameters.has_key("ComovingCoordinates") and bool(self.parameters["ComovingCoordinates"]):
+ self.cosmological_simulation = 1
+ self.omega_lambda = self.parameters["CosmologyOmegaLambdaNow"]
+ self.omega_matter = self.parameters["CosmologyOmegaMatterNow"]
+ self.hubble_constant = self.parameters["CosmologyHubbleConstantNow"]
+ a_file = open(os.path.join(self.fullplotdir,'comoving_a'))
+ line = a_file.readline().strip()
+ a_file.close()
+ self.parameters["CosmologyCurrentRedshift"] = 1/float(line) - 1
+ self.current_redshift = self.parameters["CosmologyCurrentRedshift"]
+ else:
+ self.current_redshift = self.omega_lambda = self.omega_matter = \
+ self.hubble_constant = self.cosmological_simulation = 0.0
+
+ def _parse_fparameter_file(self):
+ """
+ Parses the fortran parameter file for Orion. Most of this will
+ be useless, but this is where it keeps mu = mass per
+ particle/m_hydrogen.
+ """
+ lines = open(self.fparameter_filename).readlines()
+ for line in lines:
+ if line.count("=") == 1:
+ param, vals = map(strip,map(rstrip,line.split("=")))
+ if vals.count("'") == 0:
+ t = map(float,[a.replace('D','e').replace('d','e') for a in vals.split()]) # all are floating point.
+ else:
+ t = vals.split()
+ if len(t) == 1:
+ self.fparameters[param] = t[0]
+ else:
+ self.fparameters[param] = t
+
+ def _parse_header_file(self):
+ """
+ Parses the BoxLib header file to get any parameters stored
+ there. Hierarchy information is read out of this file in
+ OrionHierarchy.
+
+ Currently, only Time is read here.
+ """
+ header_file = open(os.path.join(self.fullplotdir,'Header'))
+ lines = header_file.readlines()
+ header_file.close()
+ n_fields = int(lines[1])
+ self.current_time = float(lines[3+n_fields])
+
+
+
+ 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()
+ self._setup_nounits_units()
+ self.conversion_factors = defaultdict(lambda: 1.0)
+ self.time_units['1'] = 1
+ self.units['1'] = 1.0
+ self.units['unitary'] = 1.0 / (self.domain_right_edge - self.domain_left_edge).max()
+ for unit in sec_conversion.keys():
+ self.time_units[unit] = 1.0 / sec_conversion[unit]
+ for key in yt2orionFieldsDict:
+ self.conversion_factors[key] = 1.0
+
+ def _setup_nounits_units(self):
+ z = 0
+ mylog.warning("Setting 1.0 in code units to be 1.0 cm")
+ 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] / mpc_conversion["cm"]
+
+ @parallel_root_only
+ def print_key_parameters(self):
+ for a in ["current_time", "domain_dimensions", "domain_left_edge",
+ "domain_right_edge"]:
+ if not hasattr(self, a):
+ mylog.error("Missing %s in parameter file definition!", a)
+ continue
+ v = getattr(self, a)
+ mylog.info("Parameters: %-25s = %s", a, v)
+
diff -r a03fc16b5457aa3deba49a11277fbb41ffebd09e -r b2cb6f33ba486efad0f9a62f2c523900c899a65e yt/frontends/boxlib/definitions.py
--- /dev/null
+++ b/yt/frontends/boxlib/definitions.py
@@ -0,0 +1,77 @@
+"""
+Various definitions for various other modules and routines
+
+Author: J. S. Oishi <jsoishi at gmail.com>
+Affiliation: KIPAC/SLAC/Stanford
+Homepage: http://yt-project.org/
+License:
+ Copyright (C) 2008-2011 J.S. Oishi. All Rights Reserved.
+
+ This file is part of yt.
+
+ yt is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+"""
+from yt.funcs import *
+
+# TODO: get rid of enzo parameters we do not need
+parameterDict = {"CosmologyCurrentRedshift": float,
+ "CosmologyComovingBoxSize": float,
+ "CosmologyOmegaMatterNow": float,
+ "CosmologyOmegaLambdaNow": float,
+ "CosmologyHubbleConstantNow": float,
+ "CosmologyInitialRedshift": float,
+ "DualEnergyFormalismEta1": float,
+ "DualEnergyFormalismEta2": float,
+ "MetaDataString": str,
+ "HydroMethod": int,
+ "DualEnergyFormalism": int,
+ "InitialTime": float,
+ "ComovingCoordinates": int,
+ "DensityUnits": float,
+ "LengthUnits": float,
+ "LengthUnit": float,
+ "TemperatureUnits": float,
+ "TimeUnits": float,
+ "GravitationalConstant": float,
+ "Gamma": float,
+ "MultiSpecies": int,
+ "CompilerPrecision": str,
+ "CurrentTimeIdentifier": int,
+ "RefineBy": int,
+ "BoundaryConditionName": str,
+ "TopGridRank": int,
+ "TopGridDimensions": int,
+ "EOSSoundSpeed": float,
+ "EOSType": int,
+ "NumberOfParticleAttributes": int,
+ }
+
+
+# converts the Orion inputs file name to the Enzo/yt name expected
+# throughout the code. key is Orion name, value is Enzo/yt equivalent
+orion2enzoDict = {"amr.n_cell": "TopGridDimensions",
+ "materials.gamma": "Gamma",
+ "amr.ref_ratio": "RefineBy",
+ "castro.use_comoving": "ComovingCoordinates",
+ "castro.redshift_in": "CosmologyInitialRedshift",
+ "comoving_OmL": "CosmologyOmegaLambdaNow",
+ "comoving_OmM": "CosmologyOmegaMatterNow",
+ "comoving_h": "CosmologyHubbleConstantNow"
+ }
+
+yt2orionFieldsDict = {}
+orion2ytFieldsDict = {}
+
+orion_FAB_header_pattern = r"^FAB \(\((\d+), \([0-9 ]+\)\),\(\d+, \(([0-9 ]+)\)\)\)\(\((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\)\) (\d+)\n"
diff -r a03fc16b5457aa3deba49a11277fbb41ffebd09e -r b2cb6f33ba486efad0f9a62f2c523900c899a65e yt/frontends/boxlib/fields.py
--- /dev/null
+++ b/yt/frontends/boxlib/fields.py
@@ -0,0 +1,191 @@
+"""
+Orion-specific fields
+
+Author: J. S. Oishi <jsoishi at gmail.com>
+Affiliation: UC Berkeley
+Homepage: http://yt-project.org/
+License:
+ Copyright (C) 2008-2011 J. S. Oishi, Matthew Turk. All Rights Reserved.
+
+ This file is part of yt.
+
+ yt is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+"""
+
+import numpy as np
+
+from yt.utilities.physical_constants import \
+ mh, kboltz
+from yt.data_objects.field_info_container import \
+ FieldInfoContainer, \
+ NullFunc, \
+ TranslationFunc, \
+ FieldInfo, \
+ ValidateParameter, \
+ ValidateDataField, \
+ ValidateProperty, \
+ ValidateSpatial, \
+ ValidateGridType
+import yt.data_objects.universal_fields
+
+
+KnownOrionFields = FieldInfoContainer()
+add_orion_field = KnownOrionFields.add_field
+
+OrionFieldInfo = FieldInfoContainer.create_with_fallback(FieldInfo)
+add_field = OrionFieldInfo.add_field
+
+add_orion_field("density", function=NullFunc, take_log=True,
+ validators = [ValidateDataField("density")],
+ units=r"\rm{g}/\rm{cm}^3")
+KnownOrionFields["density"]._projected_units =r"\rm{g}/\rm{cm}^2"
+
+add_orion_field("eden", function=NullFunc, take_log=True,
+ validators = [ValidateDataField("eden")],
+ units=r"\rm{erg}/\rm{cm}^3")
+
+add_orion_field("xmom", function=NullFunc, take_log=False,
+ validators = [ValidateDataField("xmom")],
+ units=r"\rm{g}/\rm{cm^2\ s}")
+
+add_orion_field("ymom", function=NullFunc, take_log=False,
+ validators = [ValidateDataField("ymom")],
+ units=r"\rm{gm}/\rm{cm^2\ s}")
+
+add_orion_field("zmom", function=NullFunc, take_log=False,
+ validators = [ValidateDataField("zmom")],
+ units=r"\rm{g}/\rm{cm^2\ s}")
+
+translation_dict = {"x-velocity": "xvel",
+ "y-velocity": "yvel",
+ "z-velocity": "zvel",
+ "Density": "density",
+ "TotalEnergy": "eden",
+ "Temperature": "temperature",
+ "x-momentum": "xmom",
+ "y-momentum": "ymom",
+ "z-momentum": "zmom"
+ }
+
+for f,v in translation_dict.items():
+ if v not in KnownOrionFields:
+ add_orion_field(v, function=NullFunc, take_log=False,
+ validators = [ValidateDataField(v)])
+ ff = KnownOrionFields[v]
+ add_field(f, TranslationFunc(v),
+ take_log=KnownOrionFields[v].take_log,
+ units = ff._units, display_name=f)
+
+def _xVelocity(field, data):
+ """generate x-velocity from x-momentum and density
+
+ """
+ return data["xmom"]/data["density"]
+add_orion_field("x-velocity",function=_xVelocity, take_log=False,
+ units=r'\rm{cm}/\rm{s}')
+
+def _yVelocity(field,data):
+ """generate y-velocity from y-momentum and density
+
+ """
+ #try:
+ # return data["xvel"]
+ #except KeyError:
+ return data["ymom"]/data["density"]
+add_orion_field("y-velocity",function=_yVelocity, take_log=False,
+ units=r'\rm{cm}/\rm{s}')
+
+def _zVelocity(field,data):
+ """generate z-velocity from z-momentum and density
+
+ """
+ return data["zmom"]/data["density"]
+add_orion_field("z-velocity",function=_zVelocity, take_log=False,
+ units=r'\rm{cm}/\rm{s}')
+
+def _ThermalEnergy(field, data):
+ """generate thermal (gas energy). Dual Energy Formalism was
+ implemented by Stella, but this isn't how it's called, so I'll
+ leave that commented out for now.
+ """
+ #if data.pf["DualEnergyFormalism"]:
+ # return data["GasEnergy"]
+ #else:
+ return data["TotalEnergy"] - 0.5 * data["density"] * (
+ data["x-velocity"]**2.0
+ + data["y-velocity"]**2.0
+ + data["z-velocity"]**2.0 )
+add_field("ThermalEnergy", function=_ThermalEnergy,
+ units=r"\rm{ergs}/\rm{cm^3}")
+
+def _Pressure(field,data):
+ """M{(Gamma-1.0)*e, where e is thermal energy density
+ NB: this will need to be modified for radiation
+ """
+ return (data.pf["Gamma"] - 1.0)*data["ThermalEnergy"]
+add_field("Pressure", function=_Pressure, units=r"\rm{dyne}/\rm{cm}^{2}")
+
+def _Temperature(field,data):
+ return (data.pf["Gamma"]-1.0)*data.pf["mu"]*mh*data["ThermalEnergy"]/(kboltz*data["Density"])
+add_field("Temperature",function=_Temperature,units=r"\rm{Kelvin}",take_log=False)
+
+# particle fields
+
+def particle_func(p_field, dtype='float64'):
+ def _Particles(field, data):
+ io = data.hierarchy.io
+ if not data.NumberOfParticles > 0:
+ return np.array([], dtype=dtype)
+ else:
+ return io._read_particles(data, p_field).astype(dtype)
+
+ return _Particles
+
+_particle_field_list = ["mass",
+ "position_x",
+ "position_y",
+ "position_z",
+ "momentum_x",
+ "momentum_y",
+ "momentum_z",
+ "angmomen_x",
+ "angmomen_y",
+ "angmomen_z",
+ "mlast",
+ "mdeut",
+ "n",
+ "mdot",
+ "burnstate",
+ "id"]
+
+for pf in _particle_field_list:
+ pfunc = particle_func("particle_%s" % (pf))
+ add_field("particle_%s" % pf, function=pfunc,
+ validators = [ValidateSpatial(0)],
+ particle_type=True)
+
+def _ParticleMass(field, data):
+ particles = data["particle_mass"].astype('float64')
+ return particles
+
+def _ParticleMassMsun(field, data):
+ particles = data["particle_mass"].astype('float64')
+ return particles/1.989e33
+
+add_field("ParticleMass",
+ function=_ParticleMass, validators=[ValidateSpatial(0)],
+ particle_type=True)
+add_field("ParticleMassMsun",
+ function=_ParticleMassMsun, validators=[ValidateSpatial(0)],
+ particle_type=True)
diff -r a03fc16b5457aa3deba49a11277fbb41ffebd09e -r b2cb6f33ba486efad0f9a62f2c523900c899a65e yt/frontends/boxlib/io.py
--- /dev/null
+++ b/yt/frontends/boxlib/io.py
@@ -0,0 +1,160 @@
+"""
+Orion data-file handling functions
+
+Author: Matthew Turk <matthewturk at gmail.com>
+Author: J. S. Oishi <jsoishi at gmail.com>
+Affiliation: KIPAC/SLAC/Stanford
+Homepage: http://yt-project.org/
+License:
+ Copyright (C) 2007-2011 Matthew Turk. All Rights Reserved.
+
+ This file is part of yt.
+
+ yt is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+"""
+
+import os
+import numpy as np
+from yt.utilities.io_handler import \
+ BaseIOHandler
+
+from definitions import \
+ yt2orionFieldsDict
+
+class IOHandlerNative(BaseIOHandler):
+
+ _data_style = "orion_native"
+
+ def modify(self, field):
+ return field.swapaxes(0,2)
+
+ def _read_particles(self, grid, field):
+ """
+ parses the Orion Star Particle text files
+
+ """
+ index = {'particle_mass': 0,
+ 'particle_position_x': 1,
+ 'particle_position_y': 2,
+ 'particle_position_z': 3,
+ 'particle_momentum_x': 4,
+ 'particle_momentum_y': 5,
+ 'particle_momentum_z': 6,
+ 'particle_angmomen_x': 7,
+ 'particle_angmomen_y': 8,
+ 'particle_angmomen_z': 9,
+ 'particle_mlast': 10,
+ 'particle_mdeut': 11,
+ 'particle_n': 12,
+ 'particle_mdot': 13,
+ 'particle_burnstate': 14,
+ 'particle_id': 15}
+
+ def read(line, field):
+ return float(line.split(' ')[index[field]])
+
+ fn = grid.pf.fullplotdir + "/StarParticles"
+ with open(fn, 'r') as f:
+ lines = f.readlines()
+ particles = []
+ for line in lines[1:]:
+ if grid.NumberOfParticles > 0:
+ coord = read(line, "particle_position_x"), \
+ read(line, "particle_position_y"), \
+ read(line, "particle_position_z")
+ if ( (grid.LeftEdge < coord).all() and
+ (coord <= grid.RightEdge).all() ):
+ particles.append(read(line, field))
+ return np.array(particles)
+
+ def _read_data(self,grid,field):
+ """
+ reads packed multiFABs output by BoxLib in "NATIVE" format.
+
+ """
+
+ filen = os.path.expanduser(grid.filename[field])
+ off = grid._offset[field]
+ inFile = open(filen,'rb')
+ inFile.seek(off)
+ header = inFile.readline()
+ header.strip()
+
+ if grid._paranoid:
+ mylog.warn("Orion Native reader: Paranoid read mode.")
+ headerRe = re.compile(orion_FAB_header_pattern)
+ bytesPerReal,endian,start,stop,centerType,nComponents = headerRe.search(header).groups()
+
+ # we will build up a dtype string, starting with endian
+ # check endianness (this code is ugly. fix?)
+ bytesPerReal = int(bytesPerReal)
+ if bytesPerReal == int(endian[0]):
+ dtype = '<'
+ elif bytesPerReal == int(endian[-1]):
+ dtype = '>'
+ else:
+ raise ValueError("FAB header is neither big nor little endian. Perhaps the file is corrupt?")
+
+ dtype += ('f%i'% bytesPerReal) #always a floating point
+
+ # determine size of FAB
+ start = np.array(map(int,start.split(',')))
+ stop = np.array(map(int,stop.split(',')))
+
+ gridSize = stop - start + 1
+
+ error_count = 0
+ if (start != grid.start).any():
+ print "Paranoia Error: Cell_H and %s do not agree on grid start." %grid.filename
+ error_count += 1
+ if (stop != grid.stop).any():
+ print "Paranoia Error: Cell_H and %s do not agree on grid stop." %grid.filename
+ error_count += 1
+ if (gridSize != grid.ActiveDimensions).any():
+ print "Paranoia Error: Cell_H and %s do not agree on grid dimensions." %grid.filename
+ error_count += 1
+ if bytesPerReal != grid.hierarchy._bytesPerReal:
+ print "Paranoia Error: Cell_H and %s do not agree on bytes per real number." %grid.filename
+ error_count += 1
+ if (bytesPerReal == grid.hierarchy._bytesPerReal and dtype != grid.hierarchy._dtype):
+ print "Paranoia Error: Cell_H and %s do not agree on endianness." %grid.filename
+ error_count += 1
+
+ if error_count > 0:
+ raise RunTimeError("Paranoia unveiled %i differences between Cell_H and %s." % (error_count, grid.filename))
+
+ else:
+ start = grid.start_index
+ stop = grid.stop_index
+ dtype = grid.hierarchy._dtype
+ bytesPerReal = grid.hierarchy._bytesPerReal
+
+ nElements = grid.ActiveDimensions.prod()
+
+ # one field has nElements*bytesPerReal bytes and is located
+ # nElements*bytesPerReal*field_index from the offset location
+ if yt2orionFieldsDict.has_key(field):
+ fieldname = yt2orionFieldsDict[field]
+ else:
+ fieldname = field
+ field_index = grid.field_indexes[fieldname]
+ inFile.seek(int(nElements*bytesPerReal*field_index),1)
+ field = np.fromfile(inFile,count=nElements,dtype=dtype)
+ field = field.reshape(grid.ActiveDimensions, order='F')
+
+ # we can/should also check against the max and min in the header file
+
+ inFile.close()
+ return field
+
diff -r a03fc16b5457aa3deba49a11277fbb41ffebd09e -r b2cb6f33ba486efad0f9a62f2c523900c899a65e yt/frontends/boxlib/setup.py
--- /dev/null
+++ b/yt/frontends/boxlib/setup.py
@@ -0,0 +1,13 @@
+#!/usr/bin/env python
+import setuptools
+import os
+import sys
+import os.path
+
+
+def configuration(parent_package='', top_path=None):
+ from numpy.distutils.misc_util import Configuration
+ config = Configuration('orion', parent_package, top_path)
+ config.make_config_py() # installs __config__.py
+ #config.make_svn_version_py()
+ return config
https://bitbucket.org/yt_analysis/yt/commits/030017de7f82/
Changeset: 030017de7f82
Branch: yt
User: MatthewTurk
Date: 2013-05-15 19:28:19
Summary: Splitting Boxlib / Orion readers
Affected #: 1 file
diff -r b2cb6f33ba486efad0f9a62f2c523900c899a65e -r 030017de7f82cda599a2a16b00fadc9f67c34fb9 yt/frontends/boxlib/io.py
--- a/yt/frontends/boxlib/io.py
+++ b/yt/frontends/boxlib/io.py
@@ -32,52 +32,16 @@
from definitions import \
yt2orionFieldsDict
-class IOHandlerNative(BaseIOHandler):
+class IOHandlerBoxlib(BaseIOHandler):
- _data_style = "orion_native"
+ _data_style = "boxlib_base"
def modify(self, field):
return field.swapaxes(0,2)
- def _read_particles(self, grid, field):
- """
- parses the Orion Star Particle text files
-
- """
- index = {'particle_mass': 0,
- 'particle_position_x': 1,
- 'particle_position_y': 2,
- 'particle_position_z': 3,
- 'particle_momentum_x': 4,
- 'particle_momentum_y': 5,
- 'particle_momentum_z': 6,
- 'particle_angmomen_x': 7,
- 'particle_angmomen_y': 8,
- 'particle_angmomen_z': 9,
- 'particle_mlast': 10,
- 'particle_mdeut': 11,
- 'particle_n': 12,
- 'particle_mdot': 13,
- 'particle_burnstate': 14,
- 'particle_id': 15}
-
def read(line, field):
return float(line.split(' ')[index[field]])
- fn = grid.pf.fullplotdir + "/StarParticles"
- with open(fn, 'r') as f:
- lines = f.readlines()
- particles = []
- for line in lines[1:]:
- if grid.NumberOfParticles > 0:
- coord = read(line, "particle_position_x"), \
- read(line, "particle_position_y"), \
- read(line, "particle_position_z")
- if ( (grid.LeftEdge < coord).all() and
- (coord <= grid.RightEdge).all() ):
- particles.append(read(line, field))
- return np.array(particles)
-
def _read_data(self,grid,field):
"""
reads packed multiFABs output by BoxLib in "NATIVE" format.
@@ -158,3 +122,56 @@
inFile.close()
return field
+class IOHandlerOrion(IOHandlerBoxlib):
+ _data_style = "orion_native"
+
+ def _read_particles(self, grid, field):
+ """
+ parses the Orion Star Particle text files
+
+ """
+ index = {'particle_mass': 0,
+ 'particle_position_x': 1,
+ 'particle_position_y': 2,
+ 'particle_position_z': 3,
+ 'particle_momentum_x': 4,
+ 'particle_momentum_y': 5,
+ 'particle_momentum_z': 6,
+ 'particle_angmomen_x': 7,
+ 'particle_angmomen_y': 8,
+ 'particle_angmomen_z': 9,
+ 'particle_mlast': 10,
+ 'particle_mdeut': 11,
+ 'particle_n': 12,
+ 'particle_mdot': 13,
+ 'particle_burnstate': 14,
+ 'particle_id': 15}
+
+ fn = grid.pf.fullplotdir + "/StarParticles"
+ with open(fn, 'r') as f:
+ lines = f.readlines()
+ particles = []
+ for line in lines[1:]:
+ if grid.NumberOfParticles > 0:
+ coord = read(line, "particle_position_x"), \
+ read(line, "particle_position_y"), \
+ read(line, "particle_position_z")
+ if ( (grid.LeftEdge < coord).all() and
+ (coord <= grid.RightEdge).all() ):
+ particles.append(read(line, field))
+ return np.array(particles)
+
+class IOHandlerCastro(IOHandlerBoxlib):
+ _data_style = "castro_native"
+
+ def _read_particle_field(self, grid, field):
+ offset = grid._particle_offset
+ filen = os.path.expanduser(grid.particle_filename)
+ off = grid._particle_offset
+ tr = np.zeros(grid.NumberOfParticles, dtype='float64')
+ read_castro_particles(filen, off,
+ castro_particle_field_names.index(field),
+ len(castro_particle_field_names),
+ tr)
+ return tr
+
https://bitbucket.org/yt_analysis/yt/commits/0a9cf8e7251c/
Changeset: 0a9cf8e7251c
Branch: yt
User: MatthewTurk
Date: 2013-05-15 21:15:06
Summary: Rename a few things and strip out redundancies.
Orion can now be read in with this.
Affected #: 3 files
diff -r 030017de7f82cda599a2a16b00fadc9f67c34fb9 -r 0a9cf8e7251ceeaa233dce48a9c120a0b1066e7a yt/frontends/boxlib/data_structures.py
--- a/yt/frontends/boxlib/data_structures.py
+++ b/yt/frontends/boxlib/data_structures.py
@@ -1,5 +1,5 @@
"""
-Data structures for Orion.
+Data structures for Boxlib Codes
Author: J. S. Oishi <jsoishi at gmail.com>
Affiliation: KIPAC/SLAC/Stanford
@@ -46,15 +46,14 @@
from .definitions import \
orion2enzoDict, \
parameterDict, \
- yt2orionFieldsDict, \
- orion_FAB_header_pattern
+ FAB_header_pattern
from .fields import \
OrionFieldInfo, \
add_orion_field, \
KnownOrionFields
-class OrionGrid(AMRGridPatch):
+class BoxlibGrid(AMRGridPatch):
_id_offset = 0
def __init__(self, LeftEdge, RightEdge, index, level, filename, offset,
dimensions, start, stop, paranoia=False, **kwargs):
@@ -107,11 +106,11 @@
self.field_data['dx'], self.field_data['dy'], self.field_data['dz'] = self.dds
def __repr__(self):
- return "OrionGrid_%04i" % (self.id)
+ return "BoxlibGrid_%04i" % (self.id)
-class OrionHierarchy(AMRHierarchy):
- grid = OrionGrid
- def __init__(self, pf, data_style='orion_native'):
+class BoxlibHierarchy(AMRHierarchy):
+ grid = BoxlibGrid
+ def __init__(self, pf, data_style='boxlib_native'):
self.field_indexes = {}
self.parameter_file = weakref.proxy(pf)
header_filename = os.path.join(pf.fullplotdir,'Header')
@@ -122,57 +121,11 @@
self.__cache_endianness(self.levels[-1].grids[-1])
AMRHierarchy.__init__(self,pf, self.data_style)
self._populate_hierarchy()
- self._read_particles()
+ #self._read_particles()
- def _read_particles(self):
- """
- reads in particles and assigns them to grids. Will search for
- Star particles, then sink particles if no star particle file
- is found, and finally will simply note that no particles are
- found if neither works. To add a new Orion particle type,
- simply add it to the if/elif/else block.
-
- """
- self.grid_particle_count = np.zeros(len(self.grids))
-
- for particle_filename in ["StarParticles", "SinkParticles"]:
- fn = os.path.join(self.pf.fullplotdir, particle_filename)
- if os.path.exists(fn): self._read_particle_file(fn)
-
- def _read_particle_file(self, fn):
- """actually reads the orion particle data file itself.
-
- """
- if not os.path.exists(fn): return
- with open(fn, 'r') as f:
- lines = f.readlines()
- self.num_stars = int(lines[0].strip()[0])
- for line in lines[1:]:
- particle_position_x = float(line.split(' ')[1])
- particle_position_y = float(line.split(' ')[2])
- particle_position_z = float(line.split(' ')[3])
- coord = [particle_position_x, particle_position_y, particle_position_z]
- # for each particle, determine which grids contain it
- # copied from object_finding_mixin.py
- mask=np.ones(self.num_grids)
- for i in xrange(len(coord)):
- np.choose(np.greater(self.grid_left_edge[:,i],coord[i]), (mask,0), mask)
- np.choose(np.greater(self.grid_right_edge[:,i],coord[i]), (0,mask), mask)
- ind = np.where(mask == 1)
- selected_grids = self.grids[ind]
- # in orion, particles always live on the finest level.
- # so, we want to assign the particle to the finest of
- # the grids we just found
- if len(selected_grids) != 0:
- grid = sorted(selected_grids, key=lambda grid: grid.Level)[-1]
- ind = np.where(self.grids == grid)[0][0]
- self.grid_particle_count[ind] += 1
- self.grids[ind].NumberOfParticles += 1
- return True
-
def readGlobalHeader(self,filename,paranoid_read):
"""
- read the global header file for an Orion plotfile output.
+ read the global header file for an Boxlib plotfile output.
"""
counter = 0
header_file = open(filename,'r')
@@ -248,7 +201,7 @@
counter += 1
nsteps = int(self.__global_header_lines[counter])
counter += 1
- self.levels.append(OrionLevel(lev,ngrids))
+ self.levels.append(BoxlibLevel(lev,ngrids))
# open level header, extract file names and offsets for
# each grid
# read slightly out of order here: at the end of the lo,hi
@@ -316,7 +269,7 @@
"""
Cache the endianness and bytes perreal of the grids by using a
test grid and assuming that all grids have the same
- endianness. This is a pretty safe assumption since Orion uses
+ endianness. This is a pretty safe assumption since Boxlib uses
one file per processor, and if you're running on a cluster
with different endian processors, then you're on your own!
"""
@@ -326,8 +279,7 @@
inFile.close()
header.strip()
- # parse it. the patter is in OrionDefs.py
- headerRe = re.compile(orion_FAB_header_pattern)
+ headerRe = re.compile(FAB_header_pattern)
bytesPerReal,endian,start,stop,centerType,nComponents = headerRe.search(header).groups()
self._bytesPerReal = int(bytesPerReal)
if self._bytesPerReal == int(endian[0]):
@@ -389,7 +341,6 @@
dd = self._get_data_reader_dict()
dd["field_indexes"] = self.field_indexes
AMRHierarchy._setup_classes(self, dd)
- #self._add_object_class('grid', "OrionGrid", OrionGridBase, dd)
self.object_types.sort()
def _get_grid_children(self, grid):
@@ -431,19 +382,19 @@
self._data_mode = None
self._max_locations = {}
-class OrionLevel:
+class BoxlibLevel:
def __init__(self,level,ngrids):
self.level = level
self.ngrids = ngrids
self.grids = []
-class OrionStaticOutput(StaticOutput):
+class BoxlibStaticOutput(StaticOutput):
"""
This class is a stripped down class that simply reads and parses
- *filename*, without looking at the Orion hierarchy.
+ *filename*, without looking at the Boxlib hierarchy.
"""
- _hierarchy_class = OrionHierarchy
+ _hierarchy_class = BoxlibHierarchy
_fieldinfo_fallback = OrionFieldInfo
_fieldinfo_known = KnownOrionFields
@@ -621,8 +572,6 @@
self.units['unitary'] = 1.0 / (self.domain_right_edge - self.domain_left_edge).max()
for unit in sec_conversion.keys():
self.time_units[unit] = 1.0 / sec_conversion[unit]
- for key in yt2orionFieldsDict:
- self.conversion_factors[key] = 1.0
def _setup_nounits_units(self):
z = 0
@@ -643,3 +592,50 @@
v = getattr(self, a)
mylog.info("Parameters: %-25s = %s", a, v)
+class OrionHierarchy(BoxlibHierarchy):
+ def _read_particles(self):
+ """
+ reads in particles and assigns them to grids. Will search for
+ Star particles, then sink particles if no star particle file
+ is found, and finally will simply note that no particles are
+ found if neither works. To add a new Orion particle type,
+ simply add it to the if/elif/else block.
+
+ """
+ self.grid_particle_count = np.zeros(len(self.grids))
+
+ for particle_filename in ["StarParticles", "SinkParticles"]:
+ fn = os.path.join(self.pf.fullplotdir, particle_filename)
+ if os.path.exists(fn): self._read_particle_file(fn)
+
+ def _read_particle_file(self, fn):
+ """actually reads the orion particle data file itself.
+
+ """
+ if not os.path.exists(fn): return
+ with open(fn, 'r') as f:
+ lines = f.readlines()
+ self.num_stars = int(lines[0].strip()[0])
+ for line in lines[1:]:
+ particle_position_x = float(line.split(' ')[1])
+ particle_position_y = float(line.split(' ')[2])
+ particle_position_z = float(line.split(' ')[3])
+ coord = [particle_position_x, particle_position_y, particle_position_z]
+ # for each particle, determine which grids contain it
+ # copied from object_finding_mixin.py
+ mask=np.ones(self.num_grids)
+ for i in xrange(len(coord)):
+ np.choose(np.greater(self.grid_left_edge[:,i],coord[i]), (mask,0), mask)
+ np.choose(np.greater(self.grid_right_edge[:,i],coord[i]), (0,mask), mask)
+ ind = np.where(mask == 1)
+ selected_grids = self.grids[ind]
+ # in orion, particles always live on the finest level.
+ # so, we want to assign the particle to the finest of
+ # the grids we just found
+ if len(selected_grids) != 0:
+ grid = sorted(selected_grids, key=lambda grid: grid.Level)[-1]
+ ind = np.where(self.grids == grid)[0][0]
+ self.grid_particle_count[ind] += 1
+ self.grids[ind].NumberOfParticles += 1
+ return True
+
diff -r 030017de7f82cda599a2a16b00fadc9f67c34fb9 -r 0a9cf8e7251ceeaa233dce48a9c120a0b1066e7a yt/frontends/boxlib/definitions.py
--- a/yt/frontends/boxlib/definitions.py
+++ b/yt/frontends/boxlib/definitions.py
@@ -71,7 +71,4 @@
"comoving_h": "CosmologyHubbleConstantNow"
}
-yt2orionFieldsDict = {}
-orion2ytFieldsDict = {}
-
-orion_FAB_header_pattern = r"^FAB \(\((\d+), \([0-9 ]+\)\),\(\d+, \(([0-9 ]+)\)\)\)\(\((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\)\) (\d+)\n"
+FAB_header_pattern = r"^FAB \(\((\d+), \([0-9 ]+\)\),\(\d+, \(([0-9 ]+)\)\)\)\(\((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\) \((\d+,\d+,\d+)\)\) (\d+)\n"
diff -r 030017de7f82cda599a2a16b00fadc9f67c34fb9 -r 0a9cf8e7251ceeaa233dce48a9c120a0b1066e7a yt/frontends/castro/io.py
--- a/yt/frontends/castro/io.py
+++ b/yt/frontends/castro/io.py
@@ -42,17 +42,6 @@
def modify(self, field):
return field.swapaxes(0,2)
- def _read_particle_field(self, grid, field):
- offset = grid._particle_offset
- filen = os.path.expanduser(grid.particle_filename)
- off = grid._particle_offset
- tr = np.zeros(grid.NumberOfParticles, dtype='float64')
- read_castro_particles(filen, off,
- castro_particle_field_names.index(field),
- len(castro_particle_field_names),
- tr)
- return tr
-
def _read_data(self, grid, field):
"""
reads packed multiFABs output by BoxLib in "NATIVE" format.
Repository URL: https://bitbucket.org/yt_analysis/yt/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the yt-svn
mailing list