PsiData

Contents

PsiData#

PsiData(ifile, /, model='mas', **kwargs)[source]#

Open a PSI MAS or POT3D HDF file and return the appropriate data reader.

This is the recommended entry point for reading PSI model output. The function inspects the file extension and the model argument to select the correct concrete reader class, then instantiates and returns it.

Returned object

The returned object is a lazy reader that holds an open file handle. The HDF dataset is not copied into memory at construction time — data transfer happens only when read() is called (or when the underlying data handle is explicitly indexed). The object exposes:

  • read() — primary method for loading a slice of the dataset.

  • scales — Scales(r, t, p) named tuple of coordinate scale readers; each element supports the same read() interface.

  • quantity — canonical lower-case quantity string (e.g. 'br').

  • sequence — integer time-step sequence number.

  • unit — Unit for converting from code units to physical units.

  • mesh — tuple of Mesh flags describing the Yee-grid stagger position of each axis.

  • description — human-readable description of the physical quantity.

  • native_properties — full Props descriptor.

The object also supports the context-manager protocol; the file handle is closed on exit from the with block.

The read method

odata[, r, t, p] = reader.read(*args, units=None, mesh=None, scales=True)

Positional arguments — each positional argument selects elements along one spatial axis in physical (r, θ, φ) order. Supported forms:

  • Omitted / None — full axis (slice(None)).

  • int — single index; the axis is retained as a length-1 dimension.

  • slice — standard Python slice.

  • (start, stop) or (start, stop, step) — converted to a slice.

  • Ellipsis — expands to None for all remaining axes.

Keyword arguments:

  • units — output unit. String aliases: 'native' / 'code' / 'model' return values in the custom MAS code unit; 'real' / 'phys' / 'physical' decompose to CGS base units. Any other string is forwarded to Unit and must be parseable by it — this includes SI and CGS unit names ('Gauss', 'nT', 'T'), compound expressions ('km/s', 'erg/cm**2/s'), and scale prefixes ('mG', 'μT'); see astropy:unit-format for the full grammar. An Unit instance may also be passed directly.

  • mesh — target mesh stagger (MeshCodeType). Half-mesh axes that are on the main mesh in mesh are averaged to the main mesh via remesh_arr() before return.

  • scales — if True (default), return the coordinate slice for each axis as additional Quantity values after the data.

Note

PSI HDF files are written in Fortran column-major order so that numpy reads them with shape (Nφ, Nθ, Nr) — the last axis is radial. The read positional arguments and coordinate scales are always returned in physical (r, θ, φ) order regardless of storage order.

Warning

POT3D unit convention

POT3D applies no normalization to its output. The stored values are in whatever physical units the input photospheric boundary magnetogram was provided in — most commonly Gauss, but this is not encoded in the file. The default unit for POT3D quantities is therefore dimensionless_unscaled (scale factor 1), meaning that calling read(units='physical') will not perform a meaningful unit conversion unless the correct unit is supplied explicitly via the unit keyword argument at construction time:

reader = PsiData('br001.h5', model='pot3d', unit='Gauss')
data, r, t, p = reader.read()
# data.unit == u.Gauss
Parameters:
ifilePathLike

Path to the HDF4 (.hdf) or HDF5 (.h5) file.

model{‘mas’, ‘pot3d’}, optional

PSI model type. Defaults to 'mas'.

dataset_idstr, optional

Name of the dataset within the HDF file. Defaults to the PSI standard identifier for the given format. Only needed when the file uses a non-standard dataset name.

quantitystr, optional

Override the quantity name inferred from the filename or file attributes (e.g. 'br'). Must be a key in the appropriate properties mapping.

sequenceint, optional

Override the time-step sequence number inferred from the filename or file attributes.

unitstr or astropy.units.Unit, optional

Override the code-to-physical unit derived from the quantity’s Props entry. Accepts any string parseable by Unit, including named units ('Gauss', 'nT', 'km/s'), compound expressions ('erg/cm**2/s'), and scale-prefixed forms ('mG', 'μT'); see astropy:unit-format for the complete unit grammar. An Unit instance may also be passed directly.

meshMeshCodeType, optional

Override the mesh stagger inferred from the quantity’s Props entry. Useful for files whose stagger differs from the PSI convention.

Returns:
outH5MasData or H4MasData or H5Pot3dData or H4Pot3dData

An open reader object implementing the full _HdfInterface API. The concrete type depends on model and the file extension.

Raises:
ValueError

If the combination of file extension and model is not supported, or if required metadata cannot be resolved from the file, its attributes, and the supplied keyword arguments.

FileNotFoundError

If ifile does not exist on disk.

See also

astropy.units.Unit

Unit constructor — accepts strings, compound expressions, and Unit instances.

astropy.units.Quantity.to

Unit conversion used internally when a units string is supplied to read().

Examples

Read a MAS radial magnetic field file — full array with coordinate scales:

>>> from psi_io.mhd_io import PsiData
>>> reader = PsiData('br001001.h5')
>>> data, r, t, p = reader.read()
>>> data.unit                                           # code unit (MAS_b)

Convert to Gauss on the fly:

>>> data, r, t, p = reader.read(units='Gauss')

Read only the inner radial shell (indices 0–9 in r) in CGS base units, without returning coordinate arrays:

>>> data = reader.read(slice(0, 10), units='physical', scales=False)

Read a POT3D HDF4 file. The boundary magnetogram was in Gauss, so the unit must be declared explicitly — it cannot be inferred from the file:

>>> reader = PsiData('br001.hdf', model='pot3d', unit='Gauss')
>>> data, r, t, p = reader.read()

Use as a context manager to guarantee the file handle is released:

>>> with PsiData('vr001001.h5') as reader:
...     data, r, t, p = reader.read(units='km/s')

Inspect metadata without loading data:

>>> reader = PsiData('rho001001.h5')
>>> reader.quantity          # 'rho'
>>> reader.description       # 'Mass Density'
>>> reader.unit              # MAS_n (code unit)
>>> reader.mesh              # (Mesh.HALF, Mesh.HALF, Mesh.HALF)
>>> reader.shape             # (Nφ, Nθ, Nr) — numpy storage order
Contents