Merge branch 'main' of gitea.pkm.physik.tu-darmstadt.de:IPKM/python3-damaris
This commit is contained in:
2026-07-13 13:03:39 +02:00
19 changed files with 1676 additions and 98 deletions
+117 -28
View File
@@ -9,6 +9,7 @@ import numpy
import sys
import datetime
import tables
from scipy.signal import filtfilt, remez
#############################################################################
# #
# Name: Class ADC_Result #
@@ -93,8 +94,6 @@ class ADC_Result(Resultable, Drawable, DamarisFFT, Signalpath):
print("Warning ADC-Result: Tried to run \"create_data_space()\" more than once.")
return
raise ValueError("ValueError: You cant create an ADC-Result with less than 1 sample!")
for i in range(channels):
self.y.append(numpy.zeros((samples,), dtype="int16"))
@@ -147,6 +146,97 @@ class ADC_Result(Resultable, Drawable, DamarisFFT, Signalpath):
return self.is_clipped
def lowpass(self, cutoff, numtaps=None):
"""
Apply a zero-phase lowpass FIR filter to all channels.
Uses scipy.signal.filtfilt (forward-backward filtering) with a Parks-McClellan
FIR filter (remez). The forward-backward pass guarantees zero phase distortion
and doubles the effective filter order.
The stopband is automatically set to 0.1 * Nyquist beyond the cutoff. If
numtaps is not provided, it is estimated from the transition width to yield
approximately 40-50 dB stopband attenuation (typically ~40 taps).
Parameters
----------
cutoff : float
Passband edge frequency in Hz. Frequencies above this value are attenuated.
Must be a positive number below the Nyquist frequency.
numtaps : int, optional
Number of FIR taps. If omitted, computed automatically from the transition
width.
Returns
-------
ADC_Result
A new ADC_Result with the filtered data (float64 per channel).
The original object is not modified.
Raises
------
RuntimeError
If the result contains no data.
ValueError
If cutoff is invalid for the current sampling rate.
Examples
--------
>>> # Lowpass at 1 MHz (taps auto-computed, ~40)
... filtered = adc.lowpass(cutoff=1e6)
>>> # Lowpass at 500 kHz with fixed 63 taps
... filtered = adc.lowpass(cutoff=500e3, numtaps=63)
"""
if not self.contains_data():
raise RuntimeError("lowpass: no data present")
nyquist = self.sampling_rate / 2.0
cutoff = float(cutoff)
if cutoff <= 0 or cutoff >= nyquist:
raise ValueError(
f"cutoff {cutoff} Hz is outside valid range (0, {nyquist})"
)
# Transition width: 0.1 * Nyquist
trans_width = 0.1 * nyquist
stopband = min(cutoff + trans_width, nyquist * 0.99)
# Auto-estimate taps if not given
if numtaps is None:
numtaps = max(21, int(round(8.0 / (stopband - cutoff) * nyquist)))
if numtaps % 2 == 0:
numtaps += 1
coeffs = remez(numtaps, [0, cutoff, stopband, nyquist], [1, 0], fs=self.sampling_rate)
self.lock.acquire()
try:
tmp_y = []
for i in range(self.get_number_of_channels()):
data = numpy.asarray(self.y[i], dtype='float64')
filtered = filtfilt(
coeffs, [1.0], data,
padtype='odd',
padlen=(numtaps - 1) * 3
)
tmp_y.append(filtered)
r = ADC_Result(
x=self.x[:],
y=tmp_y,
index=self.index[:],
sampl_freq=self.sampling_rate,
desc=self.description,
job_id=self.job_id,
job_date=self.job_date,
)
finally:
self.lock.release()
return r
def add_sample_space(self, samples):
"Adds space for n samples, where n can also be negative (deletes space). New space is filled up with \"0\""
@@ -220,26 +310,25 @@ class ADC_Result(Resultable, Drawable, DamarisFFT, Signalpath):
destination can be a file or a filename
suitable for further processing
"""
# write sorted
the_destination=destination
if type(destination) in (str,):
the_destination=open(destination, "w")
if isinstance(destination, str):
with open(destination, "w") as f:
self._write_to_csv_dest(f, delimiter)
else:
self._write_to_csv_dest(destination, delimiter)
the_destination.write("# adc_result\n")
the_destination.write("# t y0 y1 ...\n")
def _write_to_csv_dest(self, dest, delimiter):
dest.write("# adc_result\n")
dest.write("# t y0 y1 ...\n")
self.lock.acquire()
try:
xdata=self.get_xdata()
ch_no=self.get_number_of_channels()
ydata=list(map(self.get_ydata, range(ch_no)))
#yerr=map(self.get_yerr, xrange(ch_no))
for i in range(len(xdata)):
the_destination.write("%e"%xdata[i])
dest.write("%e"%xdata[i])
for j in range(ch_no):
the_destination.write("%s%e"%(delimiter, ydata[j][i]))
the_destination.write("\n")
the_destination=None
xdata=ydata=None
dest.write("%s%e"%(delimiter, ydata[j][i]))
dest.write("\n")
finally:
self.lock.release()
@@ -249,28 +338,28 @@ class ADC_Result(Resultable, Drawable, DamarisFFT, Signalpath):
for further processing with the NMRnotebook software;
destination can be a file or a filename
"""
# write sorted
the_destination=destination
if type(destination) in (str,):
the_destination=open(destination, "w")
if isinstance(destination, str):
with open(destination, "w") as f:
self._write_to_simpson_dest(f, delimiter)
else:
self._write_to_simpson_dest(destination, delimiter)
def _write_to_simpson_dest(self, dest, delimiter):
self.lock.acquire()
try:
xdata=self.get_xdata()
the_destination.write("SIMP\n")
the_destination.write("%s%i%s"%("NP=", len(xdata), "\n"))
the_destination.write("%s%i%s"%("SW=", self.get_sampling_rate(), "\n"))
the_destination.write("TYPE=FID\n")
the_destination.write("DATA\n")
dest.write("SIMP\n")
dest.write("%s%i%s"%("NP=", len(xdata), "\n"))
dest.write("%s%i%s"%("SW=", self.get_sampling_rate(), "\n"))
dest.write("TYPE=FID\n")
dest.write("DATA\n")
ch_no=self.get_number_of_channels()
ydata=list(map(self.get_ydata, range(ch_no)))
for i in range(len(xdata)):
for j in range(ch_no):
the_destination.write("%g%s"%(ydata[j][i], delimiter))
the_destination.write("\n")
the_destination.write("END\n")
the_destination=None
xdata=ydata=None
dest.write("%g%s"%(ydata[j][i], delimiter))
dest.write("\n")
dest.write("END\n")
finally:
self.lock.release()
+29 -24
View File
@@ -248,38 +248,41 @@ class Accumulation(Errorable, Drawable, DamarisFFT, Signalpath):
writes the data to a file.
destination can be a filehandle or a filename, default sys.stdout
"""
the_destination=destination
if isinstance(destination, str):
the_destination=open(destination, "w")
with open(destination, "w") as f:
self._write_to_csv_dest(f, delimiter)
else:
self._write_to_csv_dest(destination, delimiter)
the_destination.write("# accumulation %d\n"%self.n)
def _write_to_csv_dest(self, dest, delimiter):
dest.write("# accumulation %d\n"%self.n)
self.lock.acquire()
try:
if self.common_descriptions is not None:
for (key,value) in self.common_descriptions.items():
the_destination.write("# %s : %s\n"%(key, str(value)))
the_destination.write("# t")
dest.write("# %s : %s\n"%(key, str(value)))
dest.write("# t")
ch_no=self.get_number_of_channels()
if self.use_error:
for i in range(ch_no):
the_destination.write(" ch%d_mean ch%d_err"%(i,i))
dest.write(" ch%d_mean ch%d_err"%(i,i))
else:
for i in range(ch_no):
the_destination.write(" ch%d_mean"%i)
the_destination.write("\n")
dest.write(" ch%d_mean"%i)
dest.write("\n")
xdata=self.get_xdata()
ydata=list(map(self.get_ydata, range(ch_no)))
yerr=None
if self.use_error:
yerr=list(map(self.get_yerr, range(ch_no)))
for i in range(len(xdata)):
the_destination.write("%e"%xdata[i])
dest.write("%e"%xdata[i])
for j in range(ch_no):
if self.use_error:
the_destination.write("%s%e%s%e"%(delimiter, ydata[j][i], delimiter, yerr[j][i]))
dest.write("%s%e%s%e"%(delimiter, ydata[j][i], delimiter, yerr[j][i]))
else:
the_destination.write("%s%e"%(delimiter,ydata[j][i]))
the_destination.write("\n")
dest.write("%s%e"%(delimiter,ydata[j][i]))
dest.write("\n")
finally:
self.lock.release()
@@ -290,27 +293,29 @@ class Accumulation(Errorable, Drawable, DamarisFFT, Signalpath):
for further processing with the NMRnotebook software;
destination can be a file or a filename
"""
the_destination=destination
if isinstance(destination, str):
the_destination=open(destination, "w")
with open(destination, "w") as f:
self._write_to_simpson_dest(f, delimiter, frequency)
else:
self._write_to_simpson_dest(destination, delimiter, frequency)
def _write_to_simpson_dest(self, dest, delimiter, frequency):
self.lock.acquire()
try:
xdata=self.get_xdata()
the_destination.write("SIMP\n")
the_destination.write("%s%i%s"%("NP=", len(xdata), "\n"))
the_destination.write("%s%i%s"%("SW=", self.get_sampling_rate(), "\n"))
the_destination.write("%s%i%s"%("REF=", frequency, "\n"))
the_destination.write("TYPE=FID\n")
the_destination.write("DATA\n")
dest.write("SIMP\n")
dest.write("%s%i%s"%("NP=", len(xdata), "\n"))
dest.write("%s%i%s"%("SW=", self.get_sampling_rate(), "\n"))
dest.write("%s%i%s"%("REF=", frequency, "\n"))
dest.write("TYPE=FID\n")
dest.write("DATA\n")
ch_no=self.get_number_of_channels()
ydata=list(map(self.get_ydata, range(ch_no)))
for i in range(len(xdata)):
for j in range(ch_no):
the_destination.write("%g%s"%(ydata[j][i], delimiter))
the_destination.write("\n")
the_destination.write("END\n")
the_destination.close()
dest.write("%g%s"%(ydata[j][i], delimiter))
dest.write("\n")
dest.write("END\n")
finally:
self.lock.release()
+7 -1
View File
@@ -365,7 +365,13 @@ class DataPool(collections.abc.MutableMapping):
obj = read_from_hdf(node)
if obj is not None:
self[full_path] = obj
elif damaris_type == "TemperatureResult":
from damaris.data.Temperature import read_from_hdf
obj = read_from_hdf(node)
if obj is not None:
self[full_path] = obj
# Skip further processing of this group since we handled it
continue
+20 -22
View File
@@ -1,6 +1,5 @@
import threading
import math
import types
import sys
import tables
import numpy
@@ -45,7 +44,7 @@ class AccumulatedValue:
def __add__(self,y):
new_one=AccumulatedValue()
if (type(y) is types.InstanceType and isinstance(y, AccumulatedValue)):
if isinstance(y, AccumulatedValue):
new_one.y=self.y+y.y
new_one.y2=self.y2+y.y2
new_one.n=self.n+y.n
@@ -56,7 +55,7 @@ class AccumulatedValue:
return new_one
def __iadd__(self,y):
if (type(y) is types.InstanceType and isinstance(y, AccumulatedValue)):
if isinstance(y, AccumulatedValue):
self.y+=y.y
self.y2+=y.y2
self.n+=y.n
@@ -203,28 +202,27 @@ class MeasurementResult(Drawable.Drawable, collections.UserDict):
destination can be a file or a filename
suitable for further processing
"""
the_destination=destination
file_opened = False
if type(destination) in (str,):
the_destination=open(destination, "w")
file_opened = True
if isinstance(destination, str):
with open(destination, "w") as f:
self._write_to_csv_dest(f, delimiter)
else:
self._write_to_csv_dest(destination, delimiter)
the_destination.write("# quantity:"+str(self.quantity_name)+"\n")
the_destination.write("# x y ysigma n\n")
def _write_to_csv_dest(self, dest, delimiter):
dest.write("# quantity:"+str(self.quantity_name)+"\n")
dest.write("# x y ysigma n\n")
for x in self.get_xdata():
y=self.data[x]
if type(y) in [float, int, int]:
the_destination.write("%e%s%e%s0%s1\n"%(x, delimiter, y, delimiter, delimiter))
y = self.data[x]
if isinstance(y, (float, int)):
dest.write("%e%s%e%s0%s1\n" % (x, delimiter, y, delimiter, delimiter))
else:
the_destination.write("%e%s%e%s%e%s%d\n"%(x,
delimiter,
y.mean(),
delimiter,
y.mean_error(),
delimiter,
y.n))
if file_opened:
the_destination.close()
dest.write("%e%s%e%s%e%s%d\n" % (x,
delimiter,
y.mean(),
delimiter,
y.mean_error(),
delimiter,
y.n))
def write_to_hdf(self, hdffile, where, name, title, complib=None, complevel=None):
+218 -8
View File
@@ -2,6 +2,9 @@
from .Resultable import Resultable
from .Drawable import Drawable
import numpy
import threading
import tables
#############################################################################
# #
@@ -15,17 +18,224 @@ from .Drawable import Drawable
class TemperatureResult(Resultable, Drawable):
"""
Specialised class of Resultable and Drawable
Contains recorded temperature data
Contains recorded temperature data.
Attributes:
x: List of timestamps in seconds.
y: List of temperature readings in Celsius.
setpoint: Optional list of setpoint temperatures in Celsius.
xlabel: Label for x-axis (default: "Time (s)").
ylabel: Label for y-axis (default: "Temperature (C)").
"""
def __init__(self, x = None, y = None, desc = None, job_id = None, job_date = None):
def __init__(self, x = None, y = None, setpoint = None, desc = None, job_id = None, job_date = None):
Resultable.__init__(self)
Drawable.__init__(self)
# Set default labels
self.xlabel = "Time (s)"
self.ylabel = "Temperature (C)"
# Initialize data lists
self.x = [] if x is None else x
self.y = [] if y is None else y
self.setpoint = [] if setpoint is None else setpoint
# Set metadata if provided
if desc is not None:
self.description = desc
if job_id is not None:
self.job_id = job_id
if job_date is not None:
self.job_date = job_date
# Listener management for GUI updates
self.__listeners = []
self.__listener_lock = threading.Lock()
def get_number_of_channels(self):
"""Returns 1 for temperature data."""
return 1
def register_listener(self, listener):
"""Register a listener to be notified when data changes."""
with self.__listener_lock:
if listener not in self.__listeners:
self.__listeners.append(listener)
def unregister_listener(self, listener):
"""Unregister a listener."""
with self.__listener_lock:
if listener in self.__listeners:
self.__listeners.remove(listener)
def get_xdata(self):
"""Returns the x data (timestamps) as numpy array."""
if isinstance(self.x, numpy.ndarray):
return self.x
return numpy.array(self.x, dtype=float) if self.x else numpy.array([], dtype=float)
def get_ydata(self, channel=0):
"""Returns the y data (temperature readings) as numpy array."""
if channel == 0:
if isinstance(self.y, numpy.ndarray):
return self.y
return numpy.array(self.y, dtype=float) if self.y else numpy.array([], dtype=float)
return numpy.array([], dtype=float)
def __len__(self):
"""Returns the number of data points."""
return len(self.y) if self.y else 0
def uses_statistics(self):
"""Returns False as temperature data doesn't use statistics."""
return False
def ready_for_drawing_error(self):
"""Returns False as temperature data doesn't have error bars."""
return False
def get_yerr(self, channel=0):
"""Returns empty error data."""
return numpy.array([], dtype=float)
def get_errorplotdata(self):
"""Returns empty error plot data."""
return [[], [], []]
def get_lineplotdata(self):
"""Returns line plot data (not used for temperature)."""
return [[], []]
def get_xmin(self):
"""Returns minimum of x."""
xdata = self.get_xdata()
return xdata.min() if len(xdata) > 0 else 0
def get_xmax(self):
"""Returns maximum of x."""
xdata = self.get_xdata()
return xdata.max() if len(xdata) > 0 else 0
def get_ymin(self):
"""Returns minimum of y."""
ydata = self.get_ydata(0)
return ydata.min() if len(ydata) > 0 else 0
def get_ymax(self):
"""Returns maximum of y."""
ydata = self.get_ydata(0)
return ydata.max() if len(ydata) > 0 else 0
def get_xminpos(self):
"""Returns smallest positive value of x."""
xdata = self.get_xdata()
mask = xdata > 0
if numpy.any(mask):
return xdata[mask].min()
return 0
def get_yminpos(self):
"""Returns smallest positive value of y."""
ydata = self.get_ydata(0)
mask = ydata > 0
if numpy.any(mask):
return ydata[mask].min()
return 0
def write_to_hdf(self, hdffile, where, name, title, complib=None, complevel=None):
"""
Write the temperature result to an HDF5 group.
if (x is None) and (y is None) and (desc is None) and (job_id is None) and (job_date is None):
pass
Stores x (timestamps), y (temperature readings), and setpoint arrays
as an HDF5 table, with metadata saved as group attributes.
"""
h5_table_format = {
"x": tables.Float64Col(),
"y": tables.Float64Col(),
"setpoint": tables.Float64Col(),
}
filter = None
if complib is not None:
if complevel is None:
complevel = 9
filter = tables.Filters(complevel=complevel, complib=complib, shuffle=1)
elif (x is not None) and (y is not None) and (desc is not None) and (job_id is not None) and (job_date is not None):
pass
with self.__listener_lock:
xdata = numpy.array(self.x, dtype=float) if self.x else numpy.array([], dtype=float)
ydata = numpy.array(self.y, dtype=float) if self.y else numpy.array([], dtype=float)
spdata = numpy.array(self.setpoint, dtype=float) if self.setpoint else numpy.array([], dtype=float)
else:
raise ValueError("Wrong usage of __init__!")
# Use the length of the longest array; shorter ones are padded with 0
n = max(len(xdata), len(ydata), 1)
if len(xdata) < n:
xdata = numpy.pad(xdata, (0, n - len(xdata)), constant_values=0.0)
if len(ydata) < n:
ydata = numpy.pad(ydata, (0, n - len(ydata)), constant_values=0.0)
if len(spdata) < n:
spdata = numpy.pad(spdata, (0, n - len(spdata)), constant_values=0.0)
temp_table = hdffile.create_table(
where=where, name=name,
description=h5_table_format,
title=title,
filters=filter,
expectedrows=n,
)
temp_table.flavor = "numpy"
temp_table.attrs.damaris_type = "TemperatureResult"
temp_table.attrs.xlabel = self.xlabel
temp_table.attrs.ylabel = self.ylabel
if hasattr(self, "job_id") and self.job_id is not None:
temp_table.attrs.job_id = self.job_id
if hasattr(self, "job_date") and self.job_date is not None:
temp_table.attrs.time = "%04d%02d%02d %02d:%02d:%02d.%03d" % (
self.job_date.year,
self.job_date.month,
self.job_date.day,
self.job_date.hour,
self.job_date.minute,
self.job_date.second,
self.job_date.microsecond // 1000,
)
if hasattr(self, "description") and self.description is not None:
for (key, value) in self.description.items():
temp_table.attrs.__setattr__("description_" + key, str(value))
row = temp_table.row
for i in range(n):
row["x"] = xdata[i]
row["y"] = ydata[i]
row["setpoint"] = spdata[i]
row.append()
temp_table.flush()
def __repr__(self):
"""String representation of the temperature result."""
return f"TemperatureResult(points={len(self)}, job_id={self.job_id})"
def read_from_hdf(hdf_node):
"""
Read a TemperatureResult object from an HDF5 table node.
Returns None if the node is not suitable.
"""
if not isinstance(hdf_node, tables.Table):
return None
if hdf_node._v_attrs.damaris_type != "TemperatureResult":
return None
tr = TemperatureResult()
tr.xlabel = getattr(hdf_node._v_attrs, "xlabel", "Time (s)")
tr.ylabel = getattr(hdf_node._v_attrs, "ylabel", "Temperature (C)")
for r in hdf_node.iterrows():
tr.x.append(r["x"])
tr.y.append(r["y"])
sp = r["setpoint"]
if sp != 0.0 or len(tr.setpoint) > 0:
tr.setpoint.append(sp)
return tr