[yt-svn] commit/yt-3.0: 3 new changesets
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Fri May 17 14:58:13 PDT 2013
3 new commits in yt-3.0:
https://bitbucket.org/yt_analysis/yt-3.0/commits/1a650e92f63f/
Changeset: 1a650e92f63f
Branch: yt-3.0
User: xarthisius
Date: 2013-05-05 13:22:58
Summary: Use physical_constants module instead of hardcoded values
Affected #: 13 files
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/analysis_modules/absorption_spectrum/absorption_line.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_line.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_line.py
@@ -24,6 +24,13 @@
"""
import numpy as np
+from yt.utilities.physical_constants import \
+ charge_proton_cgs, \
+ cm_per_km, \
+ km_per_cm, \
+ mass_electron_cgs, \
+ speed_of_light_cgs
+
def voigt(a,u):
"""
@@ -167,10 +174,10 @@
"""
## constants
- me = 1.6726231e-24 / 1836. # grams mass electron
- e = 4.8032e-10 # esu
- c = 2.99792456e5 # km/s
- ccgs = c * 1.e5 # cm/s
+ me = mass_electron_cgs # grams mass electron
+ e = charge_proton_cgs # esu
+ c = speed_of_light_cgs * km_per_cm # km/s
+ ccgs = speed_of_light_cgs # cm/s
## shift lam0 by deltav
if deltav is not None:
@@ -181,7 +188,7 @@
lam1 = lam0
## conversions
- vdop = vkms * 1.e5 # in cm/s
+ vdop = vkms * cm_per_km # in cm/s
lam0cgs = lam0 / 1.e8 # rest wavelength in cm
lam1cgs = lam1 / 1.e8 # line wavelength in cm
nu1 = ccgs / lam1cgs # line freq in Hz
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/analysis_modules/halo_finding/halo_objects.py
--- a/yt/analysis_modules/halo_finding/halo_objects.py
+++ b/yt/analysis_modules/halo_finding/halo_objects.py
@@ -45,7 +45,10 @@
from yt.utilities.performance_counters import \
yt_counters, time_function
from yt.utilities.math_utils import periodic_dist, get_rotation_matrix
-from yt.utilities.physical_constants import rho_crit_now, mass_sun_cgs
+from yt.utilities.physical_constants import \
+ rho_crit_now, \
+ mass_sun_cgs, \
+ TINY
from .hop.EnzoHop import RunHOP
from .fof.EnzoFOF import RunFOF
@@ -60,8 +63,6 @@
ParallelAnalysisInterface, \
parallel_blocking_call
-TINY = 1.e-40
-
class Halo(object):
"""
A data source that returns particle information about the members of a
@@ -1428,7 +1429,7 @@
fglob = path.join(basedir, 'halos_%d.*.bin' % n)
files = glob.glob(fglob)
halos = self._get_halos_binary(files)
- #Jc = 1.98892e33/pf['mpchcm']*1e5
+ #Jc = mass_sun_cgs/ pf['mpchcm'] * 1e5
Jc = 1.0
length = 1.0 / pf['mpchcm']
conv = dict(pos = np.array([length, length, length,
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/analysis_modules/halo_mass_function/halo_mass_function.py
--- a/yt/analysis_modules/halo_mass_function/halo_mass_function.py
+++ b/yt/analysis_modules/halo_mass_function/halo_mass_function.py
@@ -31,6 +31,11 @@
ParallelDummy, \
ParallelAnalysisInterface, \
parallel_blocking_call
+from yt.utilities.physical_constants import \
+ cm_per_mpc, \
+ mass_sun_cgs, \
+ rho_crit_now
+
class HaloMassFcn(ParallelAnalysisInterface):
"""
@@ -259,7 +264,9 @@
sigma8_unnorm = math.sqrt(self.sigma_squared_of_R(R));
sigma_normalization = self.sigma8input / sigma8_unnorm;
- rho0 = self.omega_matter0 * 2.78e+11; # in units of h^2 Msolar/Mpc^3
+ # rho0 in units of h^2 Msolar/Mpc^3
+ rho0 = self.omega_matter0 * \
+ rho_crit_now * cm_per_mpc**3 / mass_sun_cgs
# spacing in mass of our sigma calculation
dm = (float(self.log_mass_max) - self.log_mass_min)/self.num_sigma_bins;
@@ -294,7 +301,9 @@
def dndm(self):
# constants - set these before calling any functions!
- rho0 = self.omega_matter0 * 2.78e+11; # in units of h^2 Msolar/Mpc^3
+ # rho0 in units of h^2 Msolar/Mpc^3
+ rho0 = self.omega_matter0 * \
+ rho_crit_now * cm_per_mpc**3 / mass_sun_cgs
self.delta_c0 = 1.69; # critical density for turnaround (Press-Schechter)
nofmz_cum = 0.0; # keep track of cumulative number density
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/analysis_modules/halo_profiler/multi_halo_profiler.py
--- a/yt/analysis_modules/halo_profiler/multi_halo_profiler.py
+++ b/yt/analysis_modules/halo_profiler/multi_halo_profiler.py
@@ -52,6 +52,9 @@
parallel_blocking_call, \
parallel_root_only, \
parallel_objects
+from yt.utilities.physical_constants import \
+ mass_sun_cgs, \
+ rho_crit_now
from yt.visualization.fixed_resolution import \
FixedResolutionBuffer
from yt.visualization.image_writer import write_image
@@ -951,12 +954,11 @@
if 'ActualOverdensity' in profile.keys():
return
- rho_crit_now = 1.8788e-29 * self.pf.hubble_constant**2 # g cm^-3
- Msun2g = 1.989e33
- rho_crit = rho_crit_now * ((1.0 + self.pf.current_redshift)**3.0)
+ rhocritnow = rho_crit_now * self.pf.hubble_constant**2 # g cm^-3
+ rho_crit = rhocritnow * ((1.0 + self.pf.current_redshift)**3.0)
if not self.use_critical_density: rho_crit *= self.pf.omega_matter
- profile['ActualOverdensity'] = (Msun2g * profile['TotalMassMsun']) / \
+ profile['ActualOverdensity'] = (mass_sun_cgs * profile['TotalMassMsun']) / \
profile['CellVolume'] / rho_crit
def _check_for_needed_profile_fields(self):
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/analysis_modules/spectral_integrator/spectral_frequency_integrator.py
--- a/yt/analysis_modules/spectral_integrator/spectral_frequency_integrator.py
+++ b/yt/analysis_modules/spectral_integrator/spectral_frequency_integrator.py
@@ -37,6 +37,9 @@
from yt.utilities.exceptions import YTException
from yt.utilities.linear_interpolators import \
BilinearFieldInterpolator
+from yt.utilities.physical_constants import \
+ erg_per_eV, \
+ keV_per_Hz
xray_data_version = 1
@@ -101,7 +104,7 @@
np.power(10, np.concatenate([self.log_E[:-1] - 0.5 * E_diff,
[self.log_E[-1] - 0.5 * E_diff[-1],
self.log_E[-1] + 0.5 * E_diff[-1]]]))
- self.dnu = 2.41799e17 * np.diff(self.E_bins)
+ self.dnu = keV_per_Hz * np.diff(self.E_bins)
def _get_interpolator(self, data, e_min, e_max):
r"""Create an interpolator for total emissivity in a
@@ -311,7 +314,7 @@
"""
my_si = EmissivityIntegrator(filename=filename)
- energy_erg = np.power(10, my_si.log_E) * 1.60217646e-9
+ energy_erg = np.power(10, my_si.log_E) * erg_per_eV
em_0 = my_si._get_interpolator((my_si.emissivity_primordial[..., :] / energy_erg),
e_min, e_max)
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/analysis_modules/star_analysis/sfr_spectrum.py
--- a/yt/analysis_modules/star_analysis/sfr_spectrum.py
+++ b/yt/analysis_modules/star_analysis/sfr_spectrum.py
@@ -31,9 +31,13 @@
from yt.utilities.cosmology import \
Cosmology, \
EnzoCosmology
+from yt.utilities.physical_constants import \
+ sec_per_year, \
+ speed_of_light_cgs
-YEAR = 3.155693e7 # sec / year
-LIGHT = 2.997925e10 # cm / s
+
+YEAR = sec_per_year # sec / year
+LIGHT = speed_of_light_cgs # cm / s
class StarFormationRate(object):
r"""Calculates the star formation rate for a given population of
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/analysis_modules/sunrise_export/sunrise_exporter.py
--- a/yt/analysis_modules/sunrise_export/sunrise_exporter.py
+++ b/yt/analysis_modules/sunrise_export/sunrise_exporter.py
@@ -35,6 +35,9 @@
import numpy as np
from yt.funcs import *
import yt.utilities.lib as amr_utils
+from yt.utilities.physical_constants import \
+ kpc_per_cm, \
+ sec_per_year
from yt.data_objects.universal_fields import add_field
from yt.mods import *
@@ -524,7 +527,7 @@
for ax in 'xyz']).transpose()
# Velocity is cm/s, we want it to be kpc/yr
#vel *= (pf["kpc"]/pf["cm"]) / (365*24*3600.)
- vel *= 1.02268944e-14
+ vel *= kpc_per_cm * sec_per_year
if initial_mass is None:
#in solar masses
initial_mass = dd["particle_mass_initial"][idx]*pf['Msun']
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/data_objects/derived_quantities.py
--- a/yt/data_objects/derived_quantities.py
+++ b/yt/data_objects/derived_quantities.py
@@ -36,6 +36,11 @@
from yt.utilities.parallel_tools.parallel_analysis_interface import \
ParallelAnalysisInterface, parallel_objects
from yt.utilities.lib import Octree
+from yt.utilities.physical_constants import \
+ gravitational_constant_cgs,
+ mass_sun_cgs,
+ HUGE
+
__CUDA_BLOCK_SIZE = 256
@@ -263,8 +268,7 @@
M = m_enc.sum()
J = np.sqrt(((j_mag.sum(axis=0))**2.0).sum())/W
E = np.sqrt(e_term_pre.sum()/W)
- G = 6.67e-8 # cm^3 g^-1 s^-2
- spin = J * E / (M*1.989e33*G)
+ spin = J * E / (M * mass_sun_cgs * gravitational_constant_cgs)
return spin
add_quantity("BaryonSpinParameter", function=_BaryonSpinParameter,
combine_function=_combBaryonSpinParameter, n_ret=4)
@@ -348,7 +352,7 @@
# Gravitational potential energy
# We only divide once here because we have velocity in cgs, but radius is
# in code.
- G = 6.67e-8 / data.convert("cm") # cm^3 g^-1 s^-2
+ G = gravitational_constant_cgs / data.convert("cm") # cm^3 g^-1 s^-2
# Check for periodicity of the clump.
two_root = 2. * np.array(data.pf.domain_width) / np.array(data.pf.domain_dimensions)
domain_period = data.pf.domain_right_edge - data.pf.domain_left_edge
@@ -570,15 +574,15 @@
mins, maxs = [], []
for field in fields:
if data[field].size < 1:
- mins.append(1e90)
- maxs.append(-1e90)
+ mins.append(HUGE)
+ maxs.append(-HUGE)
continue
if filter is None:
if non_zero:
nz_filter = data[field]>0.0
if not nz_filter.any():
- mins.append(1e90)
- maxs.append(-1e90)
+ mins.append(HUGE)
+ maxs.append(-HUGE)
continue
else:
nz_filter = None
@@ -593,8 +597,8 @@
mins.append(np.nanmin(data[field][nz_filter]))
maxs.append(np.nanmax(data[field][nz_filter]))
else:
- mins.append(1e90)
- maxs.append(-1e90)
+ mins.append(HUGE)
+ maxs.append(-HUGE)
return len(fields), mins, maxs
def _combExtrema(data, n_fields, mins, maxs):
mins, maxs = np.atleast_2d(mins, maxs)
@@ -626,7 +630,7 @@
This function returns the location of the maximum of a set
of fields.
"""
- ma, maxi, mx, my, mz = -1e90, -1, -1, -1, -1
+ ma, maxi, mx, my, mz = -HUGE, -1, -1, -1, -1
if data[field].size > 0:
maxi = np.argmax(data[field])
ma = data[field][maxi]
@@ -644,7 +648,7 @@
This function returns the location of the minimum of a set
of fields.
"""
- ma, mini, mx, my, mz = 1e90, -1, -1, -1, -1
+ ma, mini, mx, my, mz = HUGE, -1, -1, -1, -1
if data[field].size > 0:
mini = np.argmin(data[field])
ma = data[field][mini]
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/data_objects/universal_fields.py
--- a/yt/data_objects/universal_fields.py
+++ b/yt/data_objects/universal_fields.py
@@ -48,15 +48,16 @@
NeedsParameter
from yt.utilities.physical_constants import \
- mh, \
- me, \
- sigma_thompson, \
- clight, \
- kboltz, \
- G, \
- rho_crit_now, \
- speed_of_light_cgs, \
- km_per_cm, keV_per_K
+ mass_sun_cgs, \
+ mh, \
+ me, \
+ sigma_thompson, \
+ clight, \
+ kboltz, \
+ G, \
+ rho_crit_now, \
+ speed_of_light_cgs, \
+ km_per_cm, keV_per_K
from yt.utilities.math_utils import \
get_sph_r_component, \
@@ -379,7 +380,7 @@
def _CellMass(field, data):
return data["Density"] * data["CellVolume"]
def _convertCellMassMsun(data):
- return 5.027854e-34 # g^-1
+ return 1.0 / mass_sun_cgs # g^-1
add_field("CellMass", function=_CellMass, units=r"\rm{g}")
add_field("CellMassMsun", units=r"M_{\odot}",
function=_CellMass,
@@ -458,6 +459,7 @@
# lens to source
DLS = data.pf.parameters['cosmology_calculator'].AngularDiameterDistance(
data.pf.current_redshift, data.pf.parameters['lensing_source_redshift'])
+ # TODO: convert 1.5e14 to constants
return (((DL * DLS) / DS) * (1.5e14 * data.pf.omega_matter *
(data.pf.hubble_constant / speed_of_light_cgs)**2 *
(1 + data.pf.current_redshift)))
@@ -520,7 +522,7 @@
return ((data["Density"].astype('float64')**2.0) \
*data["Temperature"]**0.5)
def _convertXRayEmissivity(data):
- return 2.168e60
+ return 2.168e60 #TODO: convert me to constants
add_field("XRayEmissivity", function=_XRayEmissivity,
convert_function=_convertXRayEmissivity,
projection_conversion="1")
@@ -927,8 +929,8 @@
add_field("MeanMolecularWeight",function=_MeanMolecularWeight,units=r"")
def _JeansMassMsun(field,data):
- MJ_constant = (((5*kboltz)/(G*mh))**(1.5)) * \
- (3/(4*3.1415926535897931))**(0.5) / 1.989e33
+ MJ_constant = (((5.0 * kboltz) / (G * mh)) ** (1.5)) * \
+ (3.0 / (4.0 * np.pi)) ** (0.5) / mass_sun_cgs
return (MJ_constant *
((data["Temperature"]/data["MeanMolecularWeight"])**(1.5)) *
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/geometry/object_finding_mixin.py
--- a/yt/geometry/object_finding_mixin.py
+++ b/yt/geometry/object_finding_mixin.py
@@ -32,6 +32,8 @@
from yt.utilities.lib import \
MatchPointsToGrids, \
GridTree
+from yt.utilities.physical_constants import \
+ HUGE
class ObjectFindingMixin(object) :
@@ -83,7 +85,7 @@
Returns (value, center) of location of minimum for a given field
"""
gI = np.where(self.grid_levels >= 0) # Slow but pedantic
- minVal = 1e100
+ minVal = HUGE
for grid in self.grids[gI[0]]:
mylog.debug("Checking %s (level %s)", grid.id, grid.Level)
val, coord = grid.find_min(field)
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/utilities/cosmology.py
--- a/yt/utilities/cosmology.py
+++ b/yt/utilities/cosmology.py
@@ -25,10 +25,15 @@
"""
import numpy as np
+from yt.utilities.physical_constants import \
+ gravitational_constant_cgs, \
+ km_per_cm, \
+ pc_per_mpc, \
+ speed_of_light_cgs
-c_kms = 2.99792458e5 # c in km/s
-G = 6.67259e-8 # cgs
-kmPerMpc = 3.08567758e19
+c_kms = speed_of_light_cgs * km_per_cm # c in km/s
+G = gravitational_constant_cgs
+kmPerMpc = km_per_pc * pc_per_mpc
class Cosmology(object):
def __init__(self, HubbleConstantNow = 71.0,
@@ -162,6 +167,7 @@
"""
# Changed 2.52e17 to 2.52e19 because H_0 is in km/s/Mpc,
# instead of 100 km/s/Mpc.
+ # TODO: Move me to physical_units
return 2.52e19 / np.sqrt(self.OmegaMatterNow) / \
self.HubbleConstantNow / np.power(1 + self.InitialRedshift,1.5)
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/utilities/physical_constants.py
--- a/yt/utilities/physical_constants.py
+++ b/yt/utilities/physical_constants.py
@@ -85,3 +85,7 @@
kboltz = boltzmann_constant_cgs
hcgs = planck_constant_cgs
sigma_thompson = cross_section_thompson_cgs
+
+# Miscellaneous
+HUGE = 1.0e90
+TINY = 1.0e-40
diff -r 346c728780cb0a4607a16608a8706178a51e1bf3 -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e yt/visualization/volume_rendering/CUDARayCast.py
--- a/yt/visualization/volume_rendering/CUDARayCast.py
+++ b/yt/visualization/volume_rendering/CUDARayCast.py
@@ -29,6 +29,8 @@
import yt.extensions.HierarchySubset as hs
import numpy as np
import h5py, time
+from yt.utilities.physical_constants import \
+ mass_hydrogen_cgs
import matplotlib;matplotlib.use("Agg");import pylab
@@ -62,7 +64,7 @@
print "Constructing transfer function."
if "Data" in fn:
- mh = np.log10(1.67e-24)
+ mh = np.log10(mass_hydrogen_cgs)
tf = ColorTransferFunction((7.5+mh, 14.0+mh))
tf.add_gaussian( 8.25+mh, 0.002, [0.2, 0.2, 0.4, 0.1])
tf.add_gaussian( 9.75+mh, 0.002, [0.0, 0.0, 0.3, 0.1])
https://bitbucket.org/yt_analysis/yt-3.0/commits/557e464b6a22/
Changeset: 557e464b6a22
Branch: yt-3.0
User: xarthisius
Date: 2013-05-08 19:32:05
Summary: Fix import
Affected #: 1 file
diff -r 1a650e92f63f47456a3dcbba4eda4cca65653c2e -r 557e464b6a22cb9a51f842775ec0f8d592312be1 yt/data_objects/derived_quantities.py
--- a/yt/data_objects/derived_quantities.py
+++ b/yt/data_objects/derived_quantities.py
@@ -37,8 +37,8 @@
ParallelAnalysisInterface, parallel_objects
from yt.utilities.lib import Octree
from yt.utilities.physical_constants import \
- gravitational_constant_cgs,
- mass_sun_cgs,
+ gravitational_constant_cgs, \
+ mass_sun_cgs, \
HUGE
https://bitbucket.org/yt_analysis/yt-3.0/commits/a3bf645f0944/
Changeset: a3bf645f0944
Branch: yt-3.0
User: xarthisius
Date: 2013-05-16 19:36:28
Summary: merging
Affected #: 6 files
diff -r 557e464b6a22cb9a51f842775ec0f8d592312be1 -r a3bf645f0944f39b2550d4c445a6d59f4a017b22 yt/frontends/artio/tests/test_outputs.py
--- /dev/null
+++ b/yt/frontends/artio/tests/test_outputs.py
@@ -0,0 +1,51 @@
+"""
+ARTIO frontend tests
+
+Author: Samuel Leitner <sam.leitner at gmail.com>
+Affiliation: University of Maryland College Park
+Homepage: http://yt-project.org/
+License:
+ Copyright (C) 2012 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 yt.testing import *
+from yt.utilities.answer_testing.framework import \
+ requires_pf, \
+ data_dir_load, \
+ PixelizedProjectionValuesTest, \
+ FieldValuesTest
+from yt.frontends.artio.api import ARTIOStaticOutput
+
+_fields = ("Temperature", "Density", "VelocityMagnitude")
+
+aiso5 = "artio/aiso_a0.9005.art"
+ at requires_pf(aiso5)
+def test_aiso5():
+ pf = data_dir_load(aiso5)
+ yield assert_equal, str(pf), "aiso_a0.9005.art"
+ dso = [ None, ("sphere", ("max", (0.1, 'unitary')))]
+ for field in _fields:
+ for axis in [0, 1, 2]:
+ for ds in dso:
+ for weight_field in [None, "Density"]:
+ yield PixelizedProjectionValuesTest(
+ aiso5, axis, field, weight_field,
+ ds)
+ yield FieldValuesTest(
+ aiso5, field, ds)
+
diff -r 557e464b6a22cb9a51f842775ec0f8d592312be1 -r a3bf645f0944f39b2550d4c445a6d59f4a017b22 yt/frontends/enzo/data_structures.py
--- a/yt/frontends/enzo/data_structures.py
+++ b/yt/frontends/enzo/data_structures.py
@@ -417,6 +417,8 @@
fields = []
for ptype in self.parameter_file["AppendActiveParticleType"]:
select_grids = self.grid_active_particle_count[ptype].flat
+ if np.any(select_grids) == False:
+ continue
gs = self.grids[select_grids > 0]
g = gs[0]
handle = h5py.File(g.filename)
diff -r 557e464b6a22cb9a51f842775ec0f8d592312be1 -r a3bf645f0944f39b2550d4c445a6d59f4a017b22 yt/frontends/gdf/data_structures.py
--- a/yt/frontends/gdf/data_structures.py
+++ b/yt/frontends/gdf/data_structures.py
@@ -39,6 +39,8 @@
StaticOutput
from yt.utilities.lib import \
get_box_grids_level
+from yt.utilities.io_handler import \
+ io_registry
from yt.utilities.definitions import \
mpc_conversion, sec_conversion
@@ -78,6 +80,10 @@
if self.pf.dimensionality < 3: self.dds[2] = 1.0
self.field_data['dx'], self.field_data['dy'], self.field_data['dz'] = self.dds
+ @property
+ def filename(self):
+ return None
+
class GDFHierarchy(GridGeometryHandler):
grid = GDFGrid
@@ -85,19 +91,19 @@
def __init__(self, pf, data_style='grid_data_format'):
self.parameter_file = weakref.proxy(pf)
self.data_style = data_style
+ self.max_level = 10 # FIXME
# for now, the hierarchy file is the parameter file!
self.hierarchy_filename = self.parameter_file.parameter_filename
self.directory = os.path.dirname(self.hierarchy_filename)
- self._fhandle = h5py.File(self.hierarchy_filename,'r')
- GridGeometryHandler.__init__(self,pf,data_style)
+ self._handle = pf._handle
+ GridGeometryHandler.__init__(self, pf, data_style)
- self._fhandle.close()
def _initialize_data_storage(self):
pass
def _detect_fields(self):
- self.field_list = self._fhandle['field_types'].keys()
+ self.field_list = self._handle['field_types'].keys()
def _setup_classes(self):
dd = self._get_data_reader_dict()
@@ -105,10 +111,10 @@
self.object_types.sort()
def _count_grids(self):
- self.num_grids = self._fhandle['/grid_parent_id'].shape[0]
+ self.num_grids = self._handle['/grid_parent_id'].shape[0]
def _parse_hierarchy(self):
- f = self._fhandle
+ f = self._handle
dxs = []
self.grids = np.empty(self.num_grids, dtype='object')
levels = (f['grid_level'][:]).copy()
@@ -139,7 +145,6 @@
for gi, g in enumerate(self.grids):
g._prepare_grid()
g._setup_dx()
-
for gi, g in enumerate(self.grids):
g.Children = self._get_grid_children(g)
for g1 in g.Children:
@@ -159,22 +164,39 @@
def _setup_derived_fields(self):
self.derived_field_list = []
+ def _get_box_grids(self, left_edge, right_edge):
+ """
+ Gets back all the grids between a left edge and right edge
+ """
+ eps = np.finfo(np.float64).eps
+ grid_i = np.where((np.all((self.grid_right_edge - left_edge) > eps, axis=1) \
+ & np.all((right_edge - self.grid_left_edge) > eps, axis=1)) == True)
+
+ return self.grids[grid_i], grid_i
+
+
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)
+ grids, grid_ind = self._get_box_grids(grid.LeftEdge, grid.RightEdge)
mask[grid_ind] = True
return [g for g in self.grids[mask] if g.Level == grid.Level + 1]
+ def _setup_data_io(self):
+ self.io = io_registry[self.data_style](self.parameter_file)
+
class GDFStaticOutput(StaticOutput):
_hierarchy_class = GDFHierarchy
_fieldinfo_fallback = GDFFieldInfo
_fieldinfo_known = KnownGDFFields
+ _handle = None
def __init__(self, filename, data_style='grid_data_format',
storage_filename = None):
- StaticOutput.__init__(self, filename, data_style)
+ if self._handle is not None: return
+ self._handle = h5py.File(filename, "r")
self.storage_filename = storage_filename
self.filename = filename
+ StaticOutput.__init__(self, filename, data_style)
def _set_units(self):
"""
@@ -202,15 +224,14 @@
except:
self.units[field_name] = 1.0
try:
- current_fields_unit = current_field.attrs['field_units'][0]
+ current_fields_unit = current_field.attrs['field_units']
except:
current_fields_unit = ""
self._fieldinfo_known.add_field(field_name, function=NullFunc, take_log=False,
units=current_fields_unit, projected_units="",
convert_function=_get_convert(field_name))
-
- self._handle.close()
- del self._handle
+ for p, v in self.units.items():
+ self.conversion_factors[p] = v
def _parse_parameter_file(self):
self._handle = h5py.File(self.parameter_filename, "r")
@@ -241,8 +262,6 @@
self.hubble_constant = self.cosmological_simulation = 0.0
self.parameters['Time'] = 1.0 # Hardcode time conversion for now.
self.parameters["HydroMethod"] = 0 # Hardcode for now until field staggering is supported.
- self._handle.close()
- del self._handle
@classmethod
def _is_valid(self, *args, **kwargs):
@@ -259,3 +278,5 @@
def __repr__(self):
return self.basename.rsplit(".", 1)[0]
+ def __del__(self):
+ self._handle.close()
diff -r 557e464b6a22cb9a51f842775ec0f8d592312be1 -r a3bf645f0944f39b2550d4c445a6d59f4a017b22 yt/frontends/gdf/io.py
--- a/yt/frontends/gdf/io.py
+++ b/yt/frontends/gdf/io.py
@@ -25,31 +25,72 @@
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.funcs import \
+ mylog
from yt.utilities.io_handler import \
- BaseIOHandler
-import h5py
+ BaseIOHandler
+
+def field_dname(grid_id, field_name):
+ return "/data/grid_%010i/%s" % (grid_id, field_name)
+
+
+# TODO all particle bits were removed
class IOHandlerGDFHDF5(BaseIOHandler):
_data_style = "grid_data_format"
_offset_string = 'data:offsets=0'
_data_string = 'data:datatype=0'
- def _field_dict(self,fhandle):
- keys = fhandle['field_types'].keys()
- val = fhandle['field_types'].keys()
- return dict(zip(keys,val))
+ def __init__(self, pf, *args, **kwargs):
+ # TODO check if _num_per_stride is needed
+ self._num_per_stride = kwargs.pop("num_per_stride", 1000000)
+ BaseIOHandler.__init__(self, *args, **kwargs)
+ self.pf = pf
+ self._handle = pf._handle
- def _read_field_names(self,grid):
- fhandle = h5py.File(grid.filename,'r')
- names = fhandle['field_types'].keys()
- fhandle.close()
- return names
- def _read_data(self,grid,field):
- fhandle = h5py.File(grid.hierarchy.hierarchy_filename,'r')
- data = (fhandle['/data/grid_%010i/'%grid.id+field][:]).copy()
- fhandle.close()
- if grid.pf.field_ordering == 1:
- return data.T
+ def _read_data_set(self, grid, field):
+ if self.pf.field_ordering == 1:
+ data = self._handle[field_dname(grid.id, field)][:].swapaxes(0, 2)
else:
- return data
+ data = self._handle[field_dname(grid.id, field)][:, :, :]
+ return data.astype("float64")
+
+ def _read_data_slice(self, grid, field, axis, coord):
+ slc = [slice(None), slice(None), slice(None)]
+ slc[axis] = slice(coord, coord + 1)
+ if self.pf.field_ordering == 1:
+ data = self._handle[field_dname(grid.id, field)][:].swapaxes(0, 2)[slc]
+ else:
+ data = self._handle[field_dname(grid.id, field)][slc]
+ return data.astype("float64")
+
+ def _read_fluid_selection(self, chunks, selector, fields, size):
+ chunks = list(chunks)
+ if any((ftype != "gas" for ftype, fname in fields)):
+ raise NotImplementedError
+ fhandle = self._handle
+ rv = {}
+ for field in fields:
+ ftype, fname = field
+ rv[field] = np.empty(
+ size, dtype=fhandle[field_dname(0, fname)].dtype)
+ ngrids = sum(len(chunk.objs) for chunk in chunks)
+ mylog.debug("Reading %s cells of %s fields in %s blocks",
+ size, [fname for ftype, fname in fields], ngrids)
+ for field in fields:
+ ftype, fname = field
+ ind = 0
+ for chunk in chunks:
+ for grid in chunk.objs:
+ mask = grid.select(selector) # caches
+ if mask is None:
+ continue
+ if self.pf.field_ordering == 1:
+ data = fhandle[field_dname(grid.id, fname)][:].swapaxes(0, 2)[mask]
+ else:
+ data = fhandle[field_dname(grid.id, fname)][mask]
+ rv[field][ind:ind + data.size] = data
+ ind += data.size
+ return rv
diff -r 557e464b6a22cb9a51f842775ec0f8d592312be1 -r a3bf645f0944f39b2550d4c445a6d59f4a017b22 yt/utilities/fortran_utils.py
--- a/yt/utilities/fortran_utils.py
+++ b/yt/utilities/fortran_utils.py
@@ -81,8 +81,6 @@
s2 = vals.pop(0)
if s1 != s2:
size = struct.calcsize(endian + "I" + "".join(n*[t]) + "I")
- print "S1 = %s ; S2 = %s ; %s %s %s = %s" % (
- s1, s2, a, n, t, size)
assert(s1 == s2)
if n == 1: v = v[0]
if type(a)==tuple:
Repository URL: https://bitbucket.org/yt_analysis/yt-3.0/
--
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