"""Main photoreflectance data class `PrsData`"""
from operator import contains
import pandas as pd
import numpy as np
import os
import matplotlib.pyplot as plt
import pickle
import datetime
import logging
log = logging.getLogger(__name__)
import scipy as scp
from typing import Callable, Optional
from devctrl.data import Data, NDArray, Columns, Transform, csv as devcsv
from devctrl.utility.file_io import sanitize_filename
from devctrl.data import spectrum
from .functions import nm_to_eV, fexp
from . import functions as f
PARTIAL_PREFIX = "PART_"
"""Prefix for partial files, containing data from a single wavelength"""
METADATA_FILENAME = PARTIAL_PREFIX + "measurement_metadata.pkl"
"""Filename for just the metadata dictionary"""
FULL_DATA_SUFFIX = "_full-data.pkl"
"""Filename suffix for the full data"""
REFLECTION = "ref"
TRANSMISSION = "tra"
ABSORPTION = "abs"
PREFER_LOAD_FROM_FULL_FILE = True
"""If True, `PrsData.load_data_from_dir` will prefer loading from a full data file, rather than combining all partial files."""
I_MIN_VOLT = 0.000
OFFSET_CALCULATION = "interpolate"
"""
Determine how the offset is calculated.
- interpolate: interpolates offset values taken at various steps in time for each wavelength data point
- otherwise, simply uses the offset measurement before the measurements for all data points
"""
[docs]
class PrsData(Data):
"""
Class managing data and metadata.
Can be initialized from data directly, or from a file or directory path.
data is a dictionary with:
- key: wavelength as int or float, or "offset_<name>"
- value: dictionary with:
- key: quantity as string: ["lock-in-R_raw", "lock-in-theta", "slock-in-theta", "lock-in-aux_raw", ... ]
- value: quantity value, or array if "_raw"
Keys for a given wavelength/offset:
- dI: delta I, the change in reflectance/transmittance/absorbance. (This is the signal amplitude "R" from the lock-in)
- I: the baseline reflectivity. (DC signal measured using Aux In of the lock-in)
- theta: phase
- theta-corr: theta between [-90, 90], assuming |theta| > 90 simply means lock-in-R is negative
- lock-in-R: magnitude value of lock-in
- lock-in-theta: phase value of lock-in
- lock-in-aux: value if lock-in aux port
- <qty>_raw: The raw measurement data (all individual samples)
- s<qty>: Standard error of the quantity
"""
def __init__(self, data_path: str|None=None,
data: dict|None=None,
metadata: dict|None=None,
data_name: str|None=None,
file_mode="r",
exp_mode=None,
):
super().__init__(data_path=data_path, file_mode=file_mode, metadata=metadata)
self.data: dict = data or {}
# INIT DIRPATH
self.dirname = None
self.dirpath = None
if self._write_mode or self._read_mode:
if not data_path:
raise ValueError("A data_path must be provided in read- or write-mode")
if not self.dirpath:
if os.path.isfile(data_path):
self.dirpath = os.path.dirname(data_path)
else:
self.dirpath = data_path
self.dirname = os.path.basename(self.dirpath)
if data is None and data_path is None:
raise ValueError("Either path or data must be defined")
if self._read_mode and data is not None:
raise ValueError("Data must be None when using read-mode")
# READ
if self._read_mode and data_path is not None: # load from file
if os.path.isdir(data_path):
self.data, md = PrsData.load_data_from_dir(data_path)
self.metadata |= md
elif os.path.isfile(data_path):
if data_path.endswith(".csv"):
self.data, md = PrsData.load_data_from_csv(data_path)
self.metadata |= md
elif data_path.endswith(".pkl") or data_path.endswith(".pkl.gz"):
self.data, md = PrsData.load_data_from_pkl(data_path)
self.metadata |= md
else:
raise NotImplementedError(f"Only .csv and .pkl files are supported")
else:
raise FileNotFoundError(f"The given path '{data_path}' is neither a file nor a directory.")
if len(self.data) == 0:
log.warning(f"Data was loaded from '{data_path}', but contains no values!")
# MODE
if exp_mode is None:
if "mode" in self.metadata:
self.mode = self.metadata["mode"]
else:
self.mode = REFLECTION
else:
self.mode = exp_mode
self.metadata["mode"] = self.mode
# NAME
self.name = ""
if data_name is None:
if "name" in self.metadata:
self.name = self.metadata["name"]
elif data_path is not None:
# set name from data path: use the basename without file ext as name
self.name = os.path.basename(data_path.removesuffix("/")).split(".")[0]
log.debug(f"Setting name from data_path to '{self.name}'")
else:
self.name = data_name
# self.wavelengths = np.empty((0,), dtype=float)
self.wavelengths = []
self.update_wavelengths()
# INIT WRITE MODE
if self._write_mode:
self._assert_directory_exists()
# calculation variables
self.I_offset = None
self._offset_values = None
self._offset_value_keys = None
self._offset_times = None
def __setitem__(self, key, value):
self.data[key] = value
if type(key) == int or type(key) == float:
# self.wavelengths = np.append(self.wavelengths, key)
self.wavelengths.append(key)
self.wavelengths.sort()
def __getitem__(self, key):
return self.data[key]
def _check_has_wavelength(self, wl):
if wl not in self.data: raise KeyError(f"No data for wavelength '{wl}' ({type(wl)})")
def update_wavelengths(self):
self.wavelengths = []
keys = list(self.data.keys())
for key in keys:
# for some reason, the wavelengths can end up as string
if type(key) == str:
try:
if "." in key:
nkey = float(key)
else:
nkey = int(key)
self.data[nkey] = self.data[key]
del self.data[key]
key = nkey
except ValueError: # likely an offset
pass
if type(key) == np.float64: key = float(key)
if type(key) == np.int64: key = int(key)
if type(key) == int or type(key) == float:
# self.wavelengths = np.append(self.wavelengths, key)
self.wavelengths.append(key)
self.wavelengths.sort()
#
# POSTPROCESSING
#
[docs]
def get_for_wl(self, wl, key) -> float:
"""
Return the value for <key> at a given wavelength <wl> (in nanometers)
:param wl: Wavelength in nm
:param key: Key of the quantity to return.
if relative_time:
- I: the DC signal (measured by lock-in aux port and calculated from the "lock-in-aux_raw" values)
- dI: sqrt(dI-X^2 + dI-Y^2)
- dI-X: delta I = I * cos(theta) where I is the I value from the lock-in (peak-peak, calculated from the "lock-in-R_raw" values)
- dI-Y: delta I = I * sin(theta) where I is the I value from the lock-in (peak-peak, calculated from the "lock-in-R_raw" values)
- theta: the phase (calculated from the "lock-in-theta_raw" values)
- dI_I, dI-X_I, dI-Y_I: dI/I each
- lock-in-<key>: mean (error) value calculated from the "lock-in-<key>_raw" values
:return: float, or np.ndarray if a _raw value was requested
"""
self._check_has_wavelength(wl)
# if the key exists (because is a _raw key or it was calculated before)
if key in self.data[wl]:
return self.data[wl][key]
# circular mean calculations
if key == "lock-in-theta":
vals = self.data[wl][f"lock-in-theta_raw"]
mean = float(scp.stats.circmean(vals, low=-180, high=180))
self.data[wl][key] = mean
return mean
elif key == "slock-in-theta":
vals = self.data[wl][f"lock-in-theta_raw"]
std = float(scp.stats.circstd(vals, low=-180, high=180))
self.data[wl][key] = std
return std
# simple mean calculations
elif f"{key}_raw" in self.data[wl]:
vals = self.data[wl][f"{key}_raw"]
mean = float(np.mean(vals))
self.data[wl][key] = mean
return mean
elif key.startswith("s") and f"{key[1:]}_raw" in self.data[wl]:
vals = self.data[wl][f"{key[1:]}_raw"]
err = float(np.std(vals))
self.data[wl][key] = err
return err
# offsets: Use different calculation (obviously we must not subtract the offset when returning the offset...)
if type(wl) == str:
which = wl.removeprefix("offset_") # which offset, either begin or end
if "theta" in key:
theta, stheta = self._get_theta_offset(which)
if key.startswith("s"): return stheta
else: return theta
elif "dI" in key:
dR, sdR = self._get_dI_offset(which)
if key.startswith("s"): return sdR
else: return dR
elif "I" in key:
R, sR = self._get_I_offset(which)
if key.startswith("s"): return sR
else: return R
else:
# calculate, store and return value
calculation_funcs = {
"I": self._calculate_I_for_wl,
"dI": self._calculate_dI_for_wl,
"dI_I": self._calculate_dI_I_for_wl,
"dI-X": self._calculate_dIX_for_wl,
"dI-X_I": self._calculate_dIX_I_for_wl,
"dI-Y": self._calculate_dIY_for_wl,
"dI-Y_I": self._calculate_dIY_I_for_wl,
"theta": self._calculate_theta_for_wl,
"theta-corr": self._calculate_theta_corr_for_wl,
}
valkey = key.removeprefix("s") # value and not error
if valkey in calculation_funcs:
# value, standard error
val, sval = calculation_funcs[key](wl)
self.data[wl][valkey] = val
self.data[wl]["s" + valkey] = sval
return self.data[wl][key]
raise KeyError(f"No '{key}' data for wavelength '{wl}' ({type(wl)}). Data keys for this wavelength are {self.data[wl].keys()}")
# OFFSETS
def _remove_calculated_offsets(self):
self._offset_value_keys = None
self._offset_times = None
self._offset_values = None
def _get_time_sorted_measurement_keys(self) -> list[str|float|int]:
return sorted(self.data.keys(), key=lambda k: self.data[k]["timestamp_start"])
def _get_time_sorted_offset_keys(self) -> list[str]:
return [k for k in self._get_time_sorted_measurement_keys() if type(k) == str and "offset" in k]
def _get_time_sorted_wavelenth_keys(self) -> list[float|int]:
return [k for k in self._get_time_sorted_measurement_keys() if type(k) != str]
def _calc_offsets(self, keys=["lock-in-R", "slock-in-R", "lock-in-aux", "slock-in-aux", "lock-in-theta", "slock-in-theta"]):
"""
Initialize offset variables for later interpolation.
This function collects offset vs. time data and stores them as class members.
"""
offset_keys = self._get_time_sorted_offset_keys()
self._offset_value_keys = keys
self._offset_times = np.array([self.data[k]["timestamp_start"] for k in offset_keys])
self._offset_values = np.empty([len(offset_keys), len(self._offset_value_keys)])
for i, k in enumerate(offset_keys):
for j, vk in enumerate(self._offset_value_keys):
self._offset_values[i,j] = self.get_for_wl(wl=k, key=vk)
log.debug(f"Gathered offsets for interpolation: {offset_keys} with normalized times: {self._offset_times}")
def _calc_offset(self, wl, key) -> float:
"""
:param wl: The wavelength key to calculate the offset for.
:param key: The offset data key, for example "lock-in-R"
"""
if OFFSET_CALCULATION == "interpolate":
if self._offset_values is None:
self._calc_offsets()
# Check if there is any offset data
if self._offset_times.shape[0] == 0:
log.warning(f"{self.dirname}: No offset data - returning offset=0")
return 0.0
# Get index to offset key
try:
key_idx = self._offset_value_keys.index(key)
except ValueError as e:
# e.add_note(f"No offset interpolation was calculated for '{key}'. Avaiable are {self._offset_value_keys}")
# raise
raise Exception(f"No offset interpolation was calculated for '{key}'. Available are {self._offset_value_keys}")
# Get timestamp of measurement
try:
wl_time = self.data[wl]["timestamp_start"]
except KeyError as e:
# e.add_note(f"Can not determine timestamp for wavelenth '{wl}'. Valid are {self._offset_all_key_times.keys()}. The error was: {e}")
# raise
raise Exception(f"Can not determine timestamp for wavelenth '{wl}'. Valid are {self._offset_all_key_times.keys()}. The error was: {e}")
offset = np.interp(wl_time, self._offset_times, self._offset_values[:,key_idx])
return offset
else:
return self.get_for_wl("offset_begin", key)
# Get the offset from an offset measurement
def _get_theta_offset(self, which="begin"):
offset_name = f"offset_{which}"
ltheta, sltheta = self.get_for_wl(offset_name, "lock-in-theta"), self.get_for_wl(offset_name, "slock-in-theta")
self.data[offset_name]["theta"] = ltheta
self.data[offset_name]["stheta"] = sltheta
return ltheta, sltheta
def _get_dI_offset(self, which="begin"):
offset_name = f"offset_{which}"
dI, sdI = self.get_for_wl(offset_name, "lock-in-R"), self.get_for_wl(offset_name, "slock-in-R")
self.data[offset_name]["dI"] = dI
self.data[offset_name]["sdI"] = sdI
return dI, sdI
def _get_I_offset(self, which="begin"):
offset_name = f"offset_{which}"
try:
I, sI = self.get_for_wl(offset_name, "lock-in-aux"), self.get_for_wl(offset_name, "slock-in-aux")
self.data[offset_name]["I"] = I
self.data[offset_name]["sI"] = sI
# I offset not always recorded
except KeyError as e:
log.debug(f"[{self.dirname}]: No I offset data: '{offset_name}'")
I, sI = 0.0, 0.0
return I, sI
# CALCULATIONS
# Calcaulate desired values from the raw data
def _calculate_I_for_wl(self, wl) -> tuple[float, float]:
"""
Calculate I
If I - I_offset is negative for any wavelength, all I - I_offset values are shifted by that amount so that the smallest value is 0.
"""
# I_off, sI_off = self.get_for_wl("offset_begin", "I"), self.get_for_wl("offset_begin", "sI")
if self.I_offset is None:
# TODO calc I offset for every wavelength too?
I_off, sI_off = self._calc_offset(wl, "lock-in-aux"), self._calc_offset(wl, "slock-in-aux")
Is = np.array([self.get_for_wl(w, "lock-in-aux") - I_off for w in self.wavelengths])
I_min = np.min(Is)
if I_min <= I_MIN_VOLT:
# store for later calculations
self.I_offset = np.abs(I_MIN_VOLT - I_min) - I_off
else:
self.I_offset = - I_off
I_li = self.get_for_wl(wl, "lock-in-aux")
sI_li = self.get_for_wl(wl, "slock-in-aux")
I = I_li + self.I_offset
sI = sI_li # + self.I_offset TODO: include offset in error?
return I, sI
def _calculate_theta_for_wl(self, wl) -> tuple[float, float]:
"""
Calculate theta, which can then be used to calculate R
:param wl: wavelength key
:return: lock-in-theta - phase_offset_deg_before and mapped to [-180, 180]
"""
ltheta, sltheta = self.get_for_wl(wl, "lock-in-theta"), self.get_for_wl(wl, "slock-in-theta")
theta_off, stheta_off = self._calc_offset(wl, "lock-in-theta"), self._calc_offset(wl, "slock-in-theta")
theta = ltheta - theta_off # TODO also use phase offset after?
# map to [-180, 180]
if np.abs(theta) > 180.0:
theta = (-1.0 if theta > 0.0 else 1.0) * (360 - np.abs(theta))
# theta = theta % 180.0 * (-1)
stheta = sltheta # TODO uncertainty of phase offset
return theta, stheta
def _calculate_theta_corr_for_wl(self, wl) -> tuple[float, float]:
"""
Calculate the offset corrected phase.
The corrected phase is obtained by taking the arctan of the offset correct dI-X divided by dI-Y, which by our assumptions does not have an offset.
~~The arctan output of [-90°, 90°] is then mapped to an angle of [-180°, 180°] (angle from the positive x-axis).~~
:param wl: wavelength key
:return: corrected theta with error
"""
dIX, sdIX = self.get_for_wl(wl, "dI-X"), self.get_for_wl(wl, "sdI-X")
dIY, sdIY = self.get_for_wl(wl, "dI-Y"), self.get_for_wl(wl, "sdI-Y")
theta_corr = np.rad2deg(np.arctan(dIY/dIX))
# to map it to [-180, 180]
# if dIX < 0 and dIY < 0:
# theta_corr = -180 + theta_corr
# elif dIX < 0:
# theta_corr = 180 + theta_corr # theta_corr will be negative
# if dIX < 0:
# theta_corr *= -1
stheta_corr = np.rad2deg( 1/((dIY/dIX)**2 + 1) * np.sqrt((2 * dIY * sdIY / dIX**2)**2 + (2 * dIY**2 * sdIX / dIX**3)**2) )
# map [-180, 180] back to [-90, 90]
# if theta_corr > 90.0:
# theta_corr = -180 + theta_corr
# elif theta_corr < -90.0:
# theta_corr = 180 + theta_corr
return theta_corr, stheta_corr
def _calculate_dI_for_wl(self, wl) -> tuple[float, float]:
"""
Calculate the magnitude delta I from the lock-in's R signal and the phase theta
dI' = sqrt((delta I_X)^2 + (delta I_Y)^2)
:param wl: wavelength
:return: delta I, standard error of delta I
"""
lR, slR = self.get_for_wl(wl, "lock-in-R"), self.get_for_wl(wl, "slock-in-R")
theta, stheta = self.get_for_wl(wl, "theta"), self.get_for_wl(wl, "stheta")
dI_off, sdI_off = self._calc_offset(wl, "lock-in-R"), self._calc_offset(wl, "slock-in-R")
# TODO use offset error
return f.dI(lR, slR, theta, stheta, dI_off, sdI_off)
def _calculate_dIX_for_wl(self, wl) -> tuple[float, float]:
"""
Calculate the in-phase delta I from the lock-in's R signal and the phase theta
dI' = lock-in-I * cos(theta) - dI_offset
:param wl: wavelength
:return: delta I_X, standard error of delta I_X
"""
lR, slR = self.get_for_wl(wl, "lock-in-R"), self.get_for_wl(wl, "slock-in-R")
theta, stheta = self.get_for_wl(wl, "theta"), self.get_for_wl(wl, "stheta")
dI_off, sdI_off = self._calc_offset(wl, "lock-in-R"), self._calc_offset(wl, "slock-in-R")
# TODO use offset error
return f.dIX(lR, slR, theta, stheta, dI_off, sdI_off)
# out-of-phase
def _calculate_dIY_for_wl(self, wl) -> tuple[float, float]:
"""
Calculate the quadrature (out-of-phase) delta I from the lock-in's R signal and the phase theta
dI' = lock-in-I * sin(theta)
:param wl: wavelength
:return: delta I_Y, standard error of delta I_Y
"""
lR, slR = self.get_for_wl(wl, "lock-in-R"), self.get_for_wl(wl, "slock-in-R")
theta, stheta = self.get_for_wl(wl, "theta"), self.get_for_wl(wl, "stheta")
return f.dIY(lR, slR, theta, stheta)
# dI/I
def _calculate_dI_I_for_wl(self, wl) -> tuple[float, float]:
"""
Calculate delta I/I
:param wl: wavelength
:return: delta I/I, standard error of delta I/I
"""
dI, sdI = self.get_for_wl(wl, "dI"), self.get_for_wl(wl, "sdI")
I, sI = self.get_for_wl(wl, "I"), self.get_for_wl(wl, "sI")
return f.dI_I(dI, sdI, I, sI)
# dI/I in-phase
def _calculate_dIX_I_for_wl(self, wl) -> tuple[float, float]:
"""
Calculate delta I_X/I
:param wl: wavelength
:return: delta I_X/I, standard error of delta I_X/I
"""
dI, sdI = self.get_for_wl(wl, "dI-X"), self.get_for_wl(wl, "sdI-X")
I, sI = self.get_for_wl(wl, "I"), self.get_for_wl(wl, "sI")
return f.dI_I(dI, sdI, I, sI)
# dI/I out-of-phase
def _calculate_dIY_I_for_wl(self, wl) -> tuple[float, float]:
"""
Calculate delta I_Y/I
:param wl: wavelength
:return: delta I_Y/I, standard error of delta I_Y/I
"""
dI, sdI = self.get_for_wl(wl, "dI-Y"), self.get_for_wl(wl, "sdI-Y")
I, sI = self.get_for_wl(wl, "I"), self.get_for_wl(wl, "sI")
return f.dI_I(dI, sdI, I, sI)
# column keys that are calculated from raw data
CALC_KEYS = ["I", "dI", "dI-X", "dI-Y", "dI_I", "dI-X_I", "dI-Y_I", "theta", "theta-corr", "lock-in-theta", "lock-in-R", "lock-in-aux"]
OLD_KEYS = ["dI-Q", "dI_I-Q", "plot-theta", "real-theta"]
CALC_KEYS += OLD_KEYS
[docs]
def remove_calculated_values(self):
"""
Delete all calculated data values
"""
for wl in self.data.keys():
keys = list(self.data[wl].keys())
for key in keys:
if key in self.CALC_KEYS or key.removeprefix("s") in self.CALC_KEYS:
del self.data[wl][key]
# log.debug(f"remove_calculated_values: deleted key '{key}' for wl '{wl}'")
self.I_offset = None
self._remove_calculated_offsets()
# these are the default columns returned by the get_spectrum_data method
# should contain all possible keys
# default_spectrum_columns=["wl", "E", "dI_I",
# "theta", "lock-in-theta", "theta-corr", "stheta", "slock-in-theta", "stheta-corr",
# "dI", "lock-in-R", "sdI", "slock-in-R",
# "I", "lock-in-aux", "sI", "slock-in-aux"
# ]
# default csv columns
# fmt: off
default_spectrum_columns=[
"wl", "E",
"dI_I", "dI-X_I", "dI-Y_I",
"sdI_I", "sdI-X_I", "sdI-Y_I",
"dI", "dI-X", "dI-Y",
"sdI", "sdI-X", "sdI-Y",
"I", "sI",
"theta-corr", "stheta-corr",
"theta", "stheta", "lock-in-theta", "slock-in-theta",
"lock-in-R", "slock-in-R",
"lock-in-aux", "slock-in-aux"
]
spectrum_columns_abs=[
"wl", "E",
"dI_I", "dI-X_I", "dI-Y_I",
"sdI_I", "sdI-X_I", "sdI-Y_I",
"dI", "dI-X", "dI-Y",
"sdI", "sdI-X", "sdI-Y",
"I", "sI",
"theta-corr", "stheta-corr",
]
spectrum_columns_reference=[
"wl", "E", "I", "sI", "lock-in-aux", "slock-in-aux"
]
# fmt: on
[docs]
def get_spectrum_data(self, wavelengths=None, columns=None, Es=None, transforms: Optional[list[Transform]]=None) -> tuple[NDArray, Columns]:
"""
Return the spectral data for the specified keys and wavelengths as numpy array
:param wavelengths: List of wavelengths, range of wavelenthgs (as tuple, with at least (min, max)), or None to use all wavelengths. Conflicts with `Es`.
:param Es: List of energies (eV), range of energies (as tuple, with at least (min, max)), or None to use all wavelengths. Conflicts with `wavelengths`
:param columns: List of column keys, or None to use default_spectrum_columns.
:param transforms: List of transformations to apply to the data
:return: numpy.ndarray where the first index is the wavelength/index and the second is the keys/columns. (wavelength=0, <keys>...=1...)
"""
if columns is None: columns = self.default_spectrum_columns
if wavelengths is not None or Es is not None:
# handle offset data requested
if wavelengths is not None and any([isinstance(wl, str) for wl in wavelengths]):
pass
else:
wavelengths = spectrum.util.get_matching_or_in_range_wavelengths(self.wavelengths, wavelengths=wavelengths, Es=Es, warn_discard=True)
else:
wavelengths = self.wavelengths
data = np.empty((len(wavelengths), len(columns)), dtype=float)
# this might be slow but it shouldnt be called often
for j, wl in enumerate(wavelengths):
for i, key in enumerate(columns):
if key == "wl":
if type(wl) == str: # cant put strings in float array
data[j][i] = -1
else:
data[j][i] = wl
elif key == "E":
if type(wl) == str: # cant put strings in float array
data[j][i] = -1
else:
data[j][i] = nm_to_eV(wl)
else:
try:
val = self.get_for_wl(wl, key)
except Exception as e:
# e.add_note(f"Error getting '{key}' data for '{wl}'")
raise KeyError(f"Error getting '{key}' data for '{wl}'") from e
data[j][i] = val
if transforms:
for t in transforms:
data, columns = t(data, columns)
return data, columns
[docs]
def to_csv(self, columns=None, sep=","):
"""
Return a csv of the spectrum data as csv
:param sep: csv separator
:param columns: List of column names. If None, uses the default_spectrum_columns
:return: csv as string
"""
if columns is None:
if len(self.wavelengths) > 0:
has_cols = set(self.data[self.wavelengths[0]].keys())
# if has every necessary column to calculate the usual spectrum columns, use that
if all(map(lambda col: col in has_cols, ["lock-in-R_raw", "lock-in-aux_raw", "lock-in-theta_raw"])):
columns = self.default_spectrum_columns
# if has only aux, its likely a reference measurement
elif "lock-in-aux_raw" in has_cols:
columns = ["wl", "E", "I", "sI"]
# if it has only calculated values, it is absorption
elif all(map(lambda col: col in has_cols, self.spectrum_columns_abs)):
columns = self.spectrum_columns_abs
else:
raise ValueError(f"No columns given, and appropriate columns could not be determined")
else:
columns = ["wl", "E"] # No wavelength data, for example for pure offset measurements
sdata, columns = self.get_spectrum_data(columns=columns)
return PrsData.get_csv(sdata, self.metadata, sep=sep, columns=columns, mode=self.mode)
def save_csv_at(self, filepath, sep=",", columns=None):
log.info(f"Writing csv to {filepath}")
with open(filepath, "w") as file:
file.write(self.to_csv(sep=sep, columns=columns))
def save_csv(self, sep=",", columns=None):
self._assert_write_mode()
"""Save the csv inside the data directory"""
filepath = os.path.join(self.dirpath, self.dirname + ".csv")
self.save_csv_at(filepath, sep, columns=columns)
# FILE IO
def _assert_write_mode(self):
if not self._write_mode:
raise RuntimeError(f"Can not write data because {__class__.__name__} is not in write mode.")
def _assert_directory_exists(self):
if not os.path.isdir(self.dirpath):
log.debug(f"Making directory: {self.dirpath}")
os.makedirs(self.dirpath)
def write_partial_file(self, key):
self._assert_write_mode()
if key not in self.data:
raise KeyError(f"Invalid key '{key}'")
filename = sanitize_filename(PARTIAL_PREFIX + str(key)) + ".pkl"
self._assert_directory_exists()
filepath = os.path.join(self.dirpath, filename)
log.info(f"Writing data '{key}' to {filepath}")
with open(filepath, "wb") as file:
pickle.dump(self.data[key], file)
[docs]
def write_full_file(self, compress=False):
"""
Write the entire data to a single pickle file in the data directory.
:param compress: If True, file is compressed using gzip
"""
self._assert_write_mode()
filename = sanitize_filename(self.dirname + FULL_DATA_SUFFIX)
if compress:
filename += ".gz"
self._assert_directory_exists()
filepath = os.path.join(self.dirpath, filename)
log.info(f"Writing data to {filepath}")
if compress:
import gzip
with gzip.open(filepath, "wb") as file:
pickle.dump((self.data, self.metadata), file)
else:
with open(filepath, "wb") as file:
pickle.dump((self.data, self.metadata), file)
def write_metadata(self):
f"""
Write the metadata to the disk as '{METADATA_FILENAME}'
"""
self._assert_write_mode()
filepath = os.path.join(self.dirpath, METADATA_FILENAME)
log.debug(f"Writing metadata to {filepath}")
with open(filepath, "wb") as file:
pickle.dump(self.metadata, file)
# STATIC CONVERTER
@staticmethod
def get_csv(data: np.ndarray, metadata: dict, columns: list, mode: str, sep=","):
csv = ""
# metadata
md_keys = list(metadata.keys())
md_keys.sort()
for k in md_keys:
v = metadata[k]
if type(v) == dict:
csv += f"# {k}:\n"
keys = list(v.keys())
keys.sort()
for kk in keys:
csv += f"# {kk}:{v[kk]}\n"
elif type(v) == list:
csv += f"# {k}:\n"
for vv in v:
csv += f"# - {vv}\n"
else:
csv += f"# {k}: {v}\n"
# data header
csv += "".join(f"{PrsData.get_column_label_with_unit(colname, mode=mode)}{sep}" for colname in columns).strip(sep) + "\n"
# data
for i in range(data.shape[0]):
csv += f"{data[i, 0]}"
for j in range(1, data.shape[1]):
csv += f"{sep}{data[i,j]}"
csv += "\n"
return csv.strip("\n")
# STATIC LOADERS
# TODO
[docs]
@staticmethod
def load_data_from_csv(filepath:str, mode=None, **csv_kw) -> tuple[dict, dict]:
"""
Loads data from a single csv file.
:param filepath: Path to the csv file.
"""
sdata, columns, metadata = devcsv.from_csv(filepath, **csv_kw, yaml_load_kw=None)
return PrsData._spectrum_data_to_internal_data(sdata, columns, mode=None), metadata
@staticmethod
def _spectrum_data_to_internal_data(sdata: NDArray, columns: Columns, mode=None) -> dict:
columns = [PrsData._get_column_key("label", col, default=col, check_other_modes=True, check_other_cls=False, equality_check="startswith") for col in columns]
log.debug(f"Columns are '{columns}'")
if not "wl" in columns:
if not "E" in columns:
raise KeyError(f"The data columns must contain either 'wl' or 'E', but are {columns}")
wls = eV_to_nm(sdata[:,columns.index("E")])
else:
wls = sdata[:,columns.index("wl")]
data = {}
for i, wl in enumerate(wls):
data[wl] = {}
for j, col in enumerate(columns):
if col in ["E", "wl"]: continue
data[wl][col] = sdata[i,j]
return data
[docs]
@classmethod
def load_data_from_pkl(cls, filepath:str) -> tuple[dict, dict]:
"""
Loads data from a single pkl file.
If the file is compressed with gzip, it must end with '.gz' extension.
Parameters
----------
:param filepath Path to the file.
:return
data
2D numpy array with shape (n, 4) where n is the number of data points.
metadata
Dictionary with metadata.
"""
metadata = {}
if filepath.endswith(".gz"):
import gzip
with gzip.open(filepath, "rb") as f:
obj = pickle.load(f)
else:
with open(filepath, "rb") as f:
obj = pickle.load(f)
if isinstance(obj, tuple):
if not len(obj) == 2:
raise ValueError(f"Pickle file is a tuple with length {len(obj)}, however it must be 2: (data, metadata)")
data = obj[0]
metadata = obj[1]
if not isinstance(data, dict):
raise ValueError(f"First object in tuple is not a dictionary but {type(data)}")
elif isinstance(obj, dict):
data = obj
else:
raise ValueError(f"Pickled object must be either dict=data or (dict=data, dict=metadata), but is of type {type(obj)}")
# must be loaded by now
if not isinstance(metadata, dict):
raise ValueError(f"Metadata is not a of type dict")
return data, metadata
[docs]
@staticmethod
def load_data_from_dir(dirpath:str) -> tuple[dict, dict]:
"""
Load prs data from a directory path.
If PREFER_LOAD_FROM_FULL_FILE is True, the data is loaded from the FULL_DATA_FILENAME, if it exists.
Otherwise all data files with the PARTIAL_PREFIX are loaded and combined.
:param dirpath Path to the data directory
:return data, metadata
"""
files = os.listdir(dirpath)
if PREFER_LOAD_FROM_FULL_FILE:
for f in files:
if f.endswith(FULL_DATA_SUFFIX[1:]) or f.endswith(FULL_DATA_SUFFIX[1:] + ".gz"): # ignore - or _
return PrsData.load_data_from_pkl(os.path.join(dirpath, f))
files.sort()
data = {}
metadata = {}
for filename in files:
filepath = os.path.join(dirpath, filename)
if filename.startswith(PARTIAL_PREFIX):
log.debug(f"Loading {filename}")
# must be first
if filename == METADATA_FILENAME: # Metadata filename must also start with FLUSH_PREFIX
with open(filepath, "rb") as file:
metadata = pickle.load(file)
elif filename.endswith(".csv"):
raise NotImplementedError(f"Partial .csv files are not supported: '{filename}'")
elif filename.endswith(".pkl"):
key = filename.strip(PARTIAL_PREFIX).strip(".pkl")
with open(filepath, "rb") as file:
val = pickle.load(file)
data[key] = val
else:
raise NotImplementedError(f"Unknown file extension for file '{filepath}'")
else:
log.debug(f"Skipping unknown file: '{filepath}'")
return data, metadata
[docs]
def plot_raw_for_wl(self, wl, what=["lock-in-R", "lock-in-aux", "lock-in-theta"], fig=None, axs=None, **plot_kw):
"""
Plot raw data against index (time) for a particular wavelength.
"""
self._check_has_wavelength(wl)
if fig is None or axs is None:
fig, axs = plt.subplots(len(what), squeeze=False) # no sharex since theta might have less points
axs = axs[:,0]
if len(what) > 1: axs[-1].set_xlabel("Index")
fig.suptitle(f"Raw data for $\\lambda = {wl}$ nm")
for i, qty in enumerate(what):
axs[i].set_ylabel(PrsData.get_column_tex_label_with_unit_and_scale(qty, scale=None))
axs[i].plot(self.data[wl][f"{qty}_raw"], **(dict(color=PrsData.get_column_color(qty, default="#444")) | plot_kw))
if not fig.get_constrained_layout():
fig.tight_layout()
return fig, axs
def has_DC(self):
data = self.data[list(self.data.keys())[0]]
return "lock-in-aux_raw" in data or "I" in data
def has_AC(self):
data = self.data[list(self.data.keys())[0]]
return "lock-in-R_raw" in data or "dI" in data
def has_wavelengths(self):
return len(self.wavelengths) > 0
# fmt:off
_COLUMN_DATA = {
'*': {
'label': {
"wl": "Wavelength [nm]",
"E": "Energy [eV]",
"theta": "theta [°]",
"stheta": "sigma(theta) [°]",
"theta-corr":"theta_corr [°]",
"stheta-corr":"sigma(theta_corr) [°]",
"lock-in-R": "lock-in-R",
"slock-in-R": "sigma(lock-in-R)",
},
'label-tex': {
"wl": r"$\lambda$",
"E": r"$E$",
"theta": r"$\theta$",
"stheta": r"$\sigma(\theta)$",
"theta-corr": r"$\theta_\text{corr}$",
"stheta-corr": r"$\sigma(\theta_\text{corr})$",
"L": r"$L$",
"dI_I": r"$\Delta I/I$",
"dI-X_I": r"$\Delta I_X/I$",
"dI-Y_I": r"$\Delta I_Y/I$",
},
'unit-tex': {
"wl": r"nm",
"E": r"eV",
"dI": r"V",
"dI-X": r"V",
"dI-Y": r"V",
"I": r"V",
"theta": r"°",
"theta-corr": r"°",
"sdI": r"V",
"sdI-X": r"V",
"sdI-Y": r"V",
"sI": r"V",
"stheta": r"°",
"stheta-corr": r"°",
"lock-in-R": "V",
"lock-in-aux": "V",
"lock-in-theta": "°",
},
'scale': {
"dI_I": 1e-5,
"dI-X_I": 1e-5,
"dI-Y_I": 1e-5,
"L": 1e-5,
"dI": 1e-6,
"dI-X": 1e-6,
"dI-Y": 1e-6,
"lock-in-R": 1e-6,
},
'color': {
key: value
for p in ['', 's']
for key, value in {
p + "dI_I": "#e31a1c",
p + "dI": "#33a02c",
p + "dI-X_I": "#c34a0c",
p + "dI-X": "#43802c",
p + "dI-Y_I": "#c34a0c",
p + "dI-Y": "#43802c",
p + "lock-in-R": "#b2df8a",
p + "I": "#1f78b4",
p + "lock-in-aux": "#a6cee3",
p + "theta-corr": "#ff7f00",
p + "theta": "#fdbf6f",
p + "lock-in-theta": "#fdbf6f",
p + "L": "purple",
}.items()
}
},
REFLECTION: {
'label': {
"dI_I": "dR/R",
"dI-X_I": "dR_X/R",
"dI-Y_I": "dR_Y/R",
"dI": "dR [V]",
"dI-X": "dR_X [V]",
"dI-Y": "dR_Y [V]",
"I": "R [V]",
"sdI_I": "sigma(dR/R)",
"sdI-X_I": "sigma(dR_X/R)",
"sdI-Y_I": "sigma(dR_Y/R)",
"sdI": "sigma(dR) [V]",
"sdI-X": "sigma(dR_X) [V]",
"sdI-Y": "sigma(dR_Y) [V]",
"sI": "sigma(R) [V]",
},
'label-tex': {
"dI_I": r"$\Delta R/R$",
"dI-X_I": r"$\Delta R_X/R$",
"dI-Y_I": r"$\Delta R_Y/R$",
"dI": r"$\Delta R$",
"dI-X": r"$\Delta R_X$",
"dI-Y": r"$\Delta R_Y$",
"I": r"$R$",
"sdI_I": r"$\sigma(\Delta R/R$)",
"sdI-X_I": r"$\sigma(\Delta R_X/R$)",
"sdI-Y_I": r"$\sigma(\Delta R_Y/R$)",
"sdI": r"$\sigma(\Delta R$)",
"sdI-X": r"$\sigma(\Delta R_X$)",
"sdI-Y": r"$\sigma(\Delta R_Y$)",
},
},
TRANSMISSION: {
'label': {
"dI_I": "dT/T",
"dI-X_I": "dT_X/T",
"dI-Y_I": "dT_Y/T",
"dI": "dT [V]",
"dI-X": "dT_X [V]",
"dI-Y": "dT_Y [V]",
"I": "T [V]",
"sdI_I": "sigma(dT/T)",
"sdI-X_I": "sigma(dT_X/T)",
"sdI-Y_I": "sigma(dT_Y/T)",
"sdI": "sigma(dT) [V]",
"sdI-X": "sigma(dT_X) [V]",
"sdI-Y": "sigma(dT_Y) [V]",
"sI": "sigma(T) [V]",
},
'label-tex': {
"dI_I": r"$\Delta T/T$",
"dI-X_I": r"$\Delta T_X/T$",
"dI-Y_I": r"$\Delta T_Y/T$",
"dI": r"$\Delta T$",
"dI-X": r"$\Delta T_X$",
"dI-Y": r"$\Delta T_Y$",
"I": r"$T$",
"sdI_I": r"$\sigma(\Delta T/T$)",
"sdI-X_I": r"$\sigma(\Delta T_X/T$)",
"sdI-Y_I": r"$\sigma(\Delta T_Y/T$)",
"sdI": r"$\sigma(\Delta T$)",
"sdI-X": r"$\sigma(\Delta T_X$)",
"sdI-Y": r"$\sigma(\Delta T_Y$)",
},
},
ABSORPTION: {
'label': {
"dI_I": "dA/A",
"dI-X_I": "dA_X/A",
"dI-Y_I": "dA_Y/A",
"dI": "dA [V]",
"dI-X": "dA_X [V]",
"dI-Y": "dA_Y [V]",
"I": "A [V]",
"sdI_I": "sigma(dA/A)",
"sdI-X_I": "sigma(dA_X/A)",
"sdI-Y_I": "sigma(dA_Y/A)",
"sdI": "sigma(dA) [V]",
"sdI-X": "sigma(dA_X) [V]",
"sdI-Y": "sigma(dA_Y) [V]",
"sI": "sigma(A) [V]",
},
'label-tex': {
"dI_I": r"$\Delta A/A$",
"dI-X_I": r"$\Delta A_X/A$",
"dI-Y_I": r"$\Delta A_Y/A$",
"dI": r"$\Delta A$",
"dI-X": r"$\Delta A_X$",
"dI-Y": r"$\Delta A_Y$",
"I": r"$A$",
"sdI_I": r"$\sigma(\Delta A/A$)",
"sdI-X_I": r"$\sigma(\Delta A_X/A$)",
"sdI-Y_I": r"$\sigma(\Delta A_Y/A$)",
"sdI": r"$\sigma(\Delta A$)",
"sdI-X": r"$\sigma(\Delta A_X$)",
"sdI-Y": r"$\sigma(\Delta A_Y$)",
},
},
# fmt: on
}