Write a NeXuS HDF5 File with links to external data¶
HDF5 files may contain links to data (or groups) in other files.
This can be used to advantage to refer to data in existing HDF5 files
and create NeXus-compliant data files. Here, we show such an example,
using the same counts
v. two_theta
data from the examples above.
We use the HDF5 external file links with NeXus data files.
f[local_addr] = h5py.ExternalLink(external_file_name, external_addr)
where f
is an open h5py.File()
object in which we will create the new link,
local_addr
is an HDF5 path address, external_file_name
is the name
(relative or absolute) of an existing HDF5 file, and external_addr
is the
HDF5 path address of the existing data in the external_file_name
to be linked.
file: external_angles.hdf5¶
Take for example, the structure of external_angles.hdf5
,
a simple HDF5 data file that contains just the two_theta
angles in an HDF5 dataset at the root level of the file.
Although this is a valid HDF5 data file, it is not a valid NeXus data file:
1angles:float64[31] = [17.926079999999999, '...', 17.92108]
2 @units = degrees
file: external_counts.hdf5¶
The data in the file external_angles.hdf5
might be referenced from
another HDF5 file (such as external_counts.hdf5
)
by an HDF5 external link. [1]
Here is an example of the structure:
1entry:NXentry
2 instrument:NXinstrument
3 detector:NXdetector
4 counts:NX_INT32[31] = [1037, '...', 1321]
5 @units = counts
6 two_theta --> file="external_angles.hdf5", path="/angles"
file: external_master.hdf5¶
A valid NeXus data file could be created that refers to the data in these files without making a copy of the data files themselves.
Note
It is necessary for all these files to be located together in the same directory for the HDF5 external file links to work properly.`
To be a valid NeXus file, it must contain a NXentry group. For the files above, it is simple to make a master file that links to the data we desire, from structure that we create. We then add the group attributes that describe the default plottable data:
data:NXdata
@signal = counts
@axes = "two_theta"
@two_theta_indices = 0
Here is (the basic structure of) external_master.hdf5
, an example:
1entry:NXentry
2@default = data
3 instrument --> file="external_counts.hdf5", path="/entry/instrument"
4 data:NXdata
5 @signal = counts
6 @axes = "two_theta"
7 @two_theta = 0
8 counts --> file="external_counts.hdf5", path="/entry/instrument/detector/counts"
9 two_theta --> file="external_angles.hdf5", path="/angles"
source code: external_example_write.py¶
Here is the complete code of a Python program, using h5py
to write a NeXus-compliant HDF5 file with links to data in other HDF5 files.
external_example_write.py: Write using HDF5 external links
1#!/usr/bin/env python
2"""
3Writes a NeXus HDF5 file using h5py with links to data in other HDF5 files.
4
5This example is based on ``writer_2_1``.
6"""
7
8from pathlib import Path
9
10import h5py
11import numpy
12
13from nexusformat.nexus import (NXdata, NXdetector, NXentry, NXfield,
14 NXinstrument, NXlink, nxopen)
15
16FILE_HDF5_MASTER = "external_master.hdf5"
17FILE_HDF5_ANGLES = "external_angles.hdf5"
18FILE_HDF5_COUNTS = "external_counts.hdf5"
19
20# ---------------------------
21
22# get some data
23filename = str(Path(__file__).absolute().parent.parent / "simple_example.dat")
24buffer = numpy.loadtxt(filename).T
25tthData = buffer[0] # float[]
26countsData = numpy.asarray(buffer[1], "int32") # int[]
27
28# put the angle data in an external (non-NeXus) HDF5 data file
29with h5py.File(FILE_HDF5_ANGLES, "w") as f:
30 ds = f.create_dataset("angles", data=tthData)
31 ds.attrs["units"] = "degrees"
32
33# put the detector counts in an external HDF5 data file
34# with *incomplete* NeXus structure (no NXdata group)
35with nxopen(FILE_HDF5_COUNTS, "w") as f:
36 f["entry"] = NXentry()
37 f["entry/instrument"] = NXinstrument()
38 f["entry/instrument/detector"] = NXdetector()
39 f["entry/instrument/detector/counts"] = NXfield(countsData, units="counts")
40 f["entry/instrument/detector/two_theta"] = NXlink("/angles",
41 FILE_HDF5_ANGLES)
42
43# create a master NeXus HDF5 file
44with nxopen(FILE_HDF5_MASTER, "w") as f:
45 f["entry"] = NXentry()
46 counts = NXlink("/entry/instrument/detector/counts", FILE_HDF5_COUNTS,
47 name="counts")
48 two_theta = NXlink("/angles", FILE_HDF5_ANGLES, name="two_theta")
49 f["entry/data"] = NXdata(counts, two_theta)
50 f["entry/data"].set_default()
51 f["entry/instrument"] = NXlink("/entry/instrument", FILE_HDF5_COUNTS)
1#!/usr/bin/env python
2"""
3Writes a NeXus HDF5 file using h5py with links to data in other HDF5 files.
4
5This example is based on ``writer_2_1``.
6"""
7
8from pathlib import Path
9import h5py
10import numpy
11
12FILE_HDF5_MASTER = "external_master.hdf5"
13FILE_HDF5_ANGLES = "external_angles.hdf5"
14FILE_HDF5_COUNTS = "external_counts.hdf5"
15
16# ---------------------------
17
18# get some data
19filename = str(Path(__file__).absolute().parent.parent / "simple_example.dat")
20buffer = numpy.loadtxt(filename).T
21tthData = buffer[0] # float[]
22countsData = numpy.asarray(buffer[1], "int32") # int[]
23
24# put the angle data in an external (non-NeXus) HDF5 data file
25with h5py.File(FILE_HDF5_ANGLES, "w") as f:
26 ds = f.create_dataset("angles", data=tthData)
27 ds.attrs["units"] = "degrees"
28
29# put the detector counts in an external HDF5 data file
30# with *incomplete* NeXus structure (no NXdata group)
31with h5py.File(FILE_HDF5_COUNTS, "w") as f:
32 nxentry = f.create_group("entry")
33 nxentry.attrs["NX_class"] = "NXentry"
34 nxinstrument = nxentry.create_group("instrument")
35 nxinstrument.attrs["NX_class"] = "NXinstrument"
36 nxdetector = nxinstrument.create_group("detector")
37 nxdetector.attrs["NX_class"] = "NXdetector"
38 ds = nxdetector.create_dataset("counts", data=countsData)
39 ds.attrs["units"] = "counts"
40 # link the "two_theta" data stored in separate file
41 local_addr = nxdetector.name + "/two_theta"
42 f[local_addr] = h5py.ExternalLink(FILE_HDF5_ANGLES, "/angles")
43
44# create a master NeXus HDF5 file
45with h5py.File(FILE_HDF5_MASTER, "w") as f:
46 f.attrs["default"] = "entry"
47 nxentry = f.create_group("entry")
48 nxentry.attrs["NX_class"] = "NXentry"
49 nxentry.attrs["default"] = "data"
50 nxdata = nxentry.create_group("data")
51 nxdata.attrs["NX_class"] = "NXdata"
52
53 # link in the signal data
54 local_addr = "/entry/data/counts"
55 external_addr = "/entry/instrument/detector/counts"
56 f[local_addr] = h5py.ExternalLink(FILE_HDF5_COUNTS, external_addr)
57 nxdata.attrs["signal"] = "counts"
58
59 # link in the axes data
60 local_addr = "/entry/data/two_theta"
61 f[local_addr] = h5py.ExternalLink(FILE_HDF5_ANGLES, "/angles")
62 nxdata.attrs["axes"] = "two_theta"
63 nxdata.attrs["two_theta_indices"] = [
64 0,
65 ]
66
67 local_addr = "/entry/instrument"
68 f[local_addr] = h5py.ExternalLink(FILE_HDF5_COUNTS, "/entry/instrument")
downloads¶
The Python code and files related to this section may be downloaded from the following table.
file |
description |
---|---|
h5dump analysis of external_angles.hdf5 |
|
HDF5 file written by external_example_write |
|
punx tree analysis of external_angles.hdf5 |
|
h5dump analysis of external_counts.hdf5 |
|
HDF5 file written by external_example_write |
|
punx tree analysis of external_counts.hdf5 |
|
h5py code to write external linking examples |
|
nexusformat code to write external linking examples |
|
h5dump analysis of external_master.hdf5 |
|
NeXus file written by external_example_write |
|
punx tree analysis of external_master.hdf5 |