Mercurial > repos > muon-spectroscopy-computational-project > larch_lcf
changeset 1:6c28339b73f7 draft
planemo upload for repository https://github.com/MaterialsGalaxy/larch-tools/tree/main/larch_lcf commit 1cf6d7160497ba58fe16a51f00d088a20934eba6
author | muon-spectroscopy-computational-project |
---|---|
date | Wed, 06 Dec 2023 13:04:01 +0000 |
parents | f59731986b61 |
children | 85284c3e1354 |
files | common.py larch_lcf.py larch_lcf.xml macros.xml test-data/PtSn_OCO_Abu_1.prj test-data/PtSn_OCO_Abu_1_29204.5.prj test-data/SnO2_29206.prj test-data/SnO2_extracted.prj test-data/Sn_foil_29200.prj test-data/Sn_foil_extracted.prj |
diffstat | 10 files changed, 170 insertions(+), 72 deletions(-) [+] |
line wrap: on
line diff
--- a/common.py Tue Nov 14 15:35:22 2023 +0000 +++ b/common.py Wed Dec 06 13:04:01 2023 +0000 @@ -7,38 +7,145 @@ def get_group(athena_group: AthenaGroup, key: str = None) -> Group: + group_keys = list(athena_group._athena_groups.keys()) if key is None: - group_keys = list(athena_group._athena_groups.keys()) key = group_keys[0] - return extract_athenagroup(athena_group._athena_groups[key]) + else: + key = key.replace("-", "_") + + try: + return extract_athenagroup(athena_group._athena_groups[key]) + except KeyError as e: + raise KeyError(f"{key} not in {group_keys}") from e + + +def read_all_groups(dat_file: str, key: str = None) -> "dict[str, Group]": + # Cannot rely on do_ABC as _larch is None + athena_group = read_athena( + dat_file, + do_preedge=False, + do_bkg=False, + do_fft=False, + ) + all_groups = {} + for key in athena_group._athena_groups.keys(): + group = get_group(athena_group, key) + pre_edge_with_defaults(group=group) + xftf_with_defaults(group=group) + all_groups[key] = group + + return all_groups + + +def read_group(dat_file: str, key: str = None): + # Cannot rely on do_ABC as _larch is None + athena_group = read_athena( + dat_file, + do_preedge=False, + do_bkg=False, + do_fft=False, + ) + group = get_group(athena_group, key) + pre_edge_with_defaults(group=group) + xftf_with_defaults(group=group) + return group -def read_group(dat_file: str, key: str = None, xftf_params: dict = None): - athena_group = read_athena(dat_file) - group = get_group(athena_group, key) - bkg_parameters = group.athena_params.bkg - print(group.athena_params.fft) - print(group.athena_params.fft.__dict__) - pre_edge( - group, - e0=bkg_parameters.e0, - pre1=bkg_parameters.pre1, - pre2=bkg_parameters.pre2, - norm1=bkg_parameters.nor1, - norm2=bkg_parameters.nor2, - nnorm=bkg_parameters.nnorm, - make_flat=bkg_parameters.flatten, +def pre_edge_with_defaults(group: Group, settings: dict = None): + merged_settings = {} + try: + bkg_parameters = group.athena_params.bkg + except AttributeError as e: + print(f"Cannot load group.athena_params.bkg from group:\n{e}") + bkg_parameters = None + + keys = ( + ("e0", "e0", None), + ("pre1", "pre1", None), + ("pre2", "pre2", None), + ("norm1", "nor1", None), + ("norm2", "nor2", None), + ("nnorm", "nnorm", None), + ("make_flat", "flatten", None), + ("step", "step", None), + # This cannot be read from file as it is not stored by Larch (0.9.71) + # ("nvict", "nvict", None), ) - autobk(group) - if xftf_params is None: - xftf(group) - else: - print(xftf_params) - xftf(group, **xftf_params) - xftf_details = Group() - setattr(xftf_details, "call_args", xftf_params) - group.xftf_details = xftf_details - return group + for key, parameters_key, default in keys: + extract_attribute( + merged_settings, key, bkg_parameters, parameters_key, default + ) + + if settings: + for k, v in settings.items(): + if k == "nvict": + print( + "WARNING: `nvict` can be used for pre-edge but is not " + "saved to file, so value used will not be accessible in " + "future operations using this Athena .prj" + ) + merged_settings[k] = v + + print(f"Pre-edge normalization with {merged_settings}") + try: + pre_edge(group, **merged_settings) + except Warning as e: + raise Warning( + "Unable to perform pre-edge fitting with:\n\n" + f"energy:\n{group.energy}\n\nmu:{group.mu}\n\n" + "Consider checking the correct columns have been extracted" + ) from e + autobk(group, pre_edge_kws=merged_settings) + + +def xftf_with_defaults(group: Group, settings: dict = None): + merged_settings = {} + try: + fft_parameters = group.athena_params.fft + except AttributeError as e: + print(f"Cannot load group.athena_params.fft from group:\n{e}") + fft_parameters = None + + keys = ( + ("kmin", "kmin", 0), + ("kmax", "kmax", 20), + ("dk", "dk", 1), + ("kweight", "kw", 2), + ("kweight", "kweight", 2), + ("window", "kwindow", "kaiser"), + ) + for key, parameters_key, default in keys: + extract_attribute( + merged_settings, key, fft_parameters, parameters_key, default + ) + + if settings: + for k, v in settings.items(): + merged_settings[k] = v + + print(f"XFTF with {merged_settings}") + xftf(group, **merged_settings) + xftf_details = Group() + setattr(xftf_details, "call_args", merged_settings) + group.xftf_details = xftf_details + + +def extract_attribute( + merged_settings: dict, + key: str, + parameters_group: Group, + parameters_key: str, + default: "str|int" = None, +): + if parameters_group is not None: + try: + merged_settings[key] = getattr(parameters_group, parameters_key) + return + except AttributeError: + pass + + if default is not None: + merged_settings[key] = default def read_groups(dat_files: "list[str]", key: str = None) -> Iterable[Group]:
--- a/larch_lcf.py Tue Nov 14 15:35:22 2023 +0000 +++ b/larch_lcf.py Wed Dec 06 13:04:01 2023 +0000 @@ -13,8 +13,8 @@ def plot( group_to_fit: Group, fit_group: Group, - energy_min: float, - energy_max: float, + x_limit_min: float, + x_limit_max: float, ): formatted_label = "" for label, weight in fit_group.weights.items(): @@ -37,7 +37,7 @@ linestyle="--", ) plt.grid(color="black", linestyle=":", linewidth=1) # show and format grid - plt.xlim(energy_min, energy_max) + plt.xlim(x_limit_min, x_limit_max) plt.xlabel("Energy (eV)") plt.ylabel("normalised x$\mu$(E)") # noqa: W605 plt.legend() @@ -67,16 +67,16 @@ set_label(component_group, component["label"]) component_groups.append(component_group) - fit_group = lincombo_fit(group_to_fit, component_groups) - print(f"Goodness of fit (rfactor): {fit_group.rfactor:.6%}") - energy_min = input_values["energy_min"] energy_max = input_values["energy_max"] - if input_values["energy_format"] == "relative": - e0 = group_to_fit.e0 - if energy_min is not None: - energy_min += e0 - if energy_max is not None: - energy_max += e0 + fit_group = lincombo_fit( + group=group_to_fit, + components=component_groups, + xmin=energy_min, + xmax=energy_max, + ) + print(f"Goodness of fit (rfactor): {fit_group.rfactor:.6%}") - plot(group_to_fit, fit_group, energy_min, energy_max) + x_limit_min = input_values["x_limit_min"] + x_limit_max = input_values["x_limit_max"] + plot(group_to_fit, fit_group, x_limit_min, x_limit_max)
--- a/larch_lcf.xml Tue Nov 14 15:35:22 2023 +0000 +++ b/larch_lcf.xml Wed Dec 06 13:04:01 2023 +0000 @@ -4,7 +4,7 @@ <!-- version of underlying tool (PEP 440) --> <token name="@TOOL_VERSION@">0.9.71</token> <!-- version of this tool wrapper (integer) --> - <token name="@WRAPPER_VERSION@">0</token> + <token name="@WRAPPER_VERSION@">1</token> <!-- citation should be updated with every underlying tool version --> <!-- typical fields to update are version, month, year, and doi --> <token name="@TOOL_CITATION@">10.1088/1742-6596/430/1/012007</token> @@ -19,6 +19,7 @@ </requirements> <required_files> <include type="literal" path="larch_lcf.py"/> + <include type="literal" path="common.py"/> </required_files> <command detect_errors="exit_code"><![CDATA[ python '${__tool_directory__}/larch_lcf.py' '$prj_file' '$inputs' @@ -33,23 +34,33 @@ <param name="label" type="text" optional="true" label="Component label" help="The label to use for this component. If unset, the label will be taken from the Athena project metadata."/> <param name="component_file" type="data" format="prj" label="Athena project to fit" help="Reference X-ray Absorption Fine Structure (XAFS) data, in Athena project format, to be used as components in the fit."/> </repeat> - <expand macro="energy_limits"/> + <!-- Energy range for LCF --> + <param name="energy_min" type="float" label="Minimum fit energy (eV)" optional="true" help="If set, only data above this value will be used for the LCF."/> + <param name="energy_max" type="float" label="Maximum fit energy (eV)" optional="true" help="If set, only data below this value will be used for the LCF."/> + <!-- Plot limits --> + <expand macro="plot_limits_energy"/> </inputs> <outputs> - <data name="plot" format="png" from_work_dir="plot.png" label="Larch LCF ${label} ${on_string}"/> + <data name="plot" format="png" from_work_dir="plot.png"/> </outputs> <tests> + <!-- 1 --> <test expect_num_outputs="1"> - <param name="prj_file" value="PtSn_OCO_Abu_1.prj"/> - <param name="component_file" value="Sn_foil_extracted.prj"/> - <param name="component_file" value="SnO2_extracted.prj"/> - <param name="energy_format" value="relative"/> - <param name="energy_max" value="10"/> + <param name="prj_file" value="PtSn_OCO_Abu_1_29204.5.prj"/> + <param name="component_file" value="Sn_foil_29200.prj"/> + <param name="component_file" value="SnO2_29206.prj"/> + <param name="energy_min" value="29190"/> + <param name="energy_max" value="29230"/> + <param name="x_limit_min" value="29190"/> + <param name="x_limit_max" value="29230"/> <output name="plot"> <assert_contents> - <has_size value="55100" delta="100"/> + <has_size value="59500" delta="100"/> </assert_contents> </output> + <assert_stdout> + <has_text text="Goodness of fit (rfactor): 0.532130%"/> + </assert_stdout> </test> </tests> <help><![CDATA[
--- a/macros.xml Tue Nov 14 15:35:22 2023 +0000 +++ b/macros.xml Wed Dec 06 13:04:01 2023 +0000 @@ -1,26 +1,6 @@ <macros> - <xml name="energy_limits"> - <param name="energy_format" type="select" display="radio" label="Energy limits" help="Whether to limit the energy relative to the absorption edge or with absolute values."> - <option value="relative" selected="true">Relative</option> - <option value="absolute">Absolute</option> - </param> - <param name="energy_min" type="float" label="Minimum energy (eV)" optional="true" help="If set, data will be cropped below this value in electron volts."/> - <param name="energy_max" type="float" label="Maximum energy (eV)" optional="true" help="If set, data will be cropped above this value in electron volts."/> - </xml> - <xml name="xftf_params"> - <param argument="kmin" type="float" value="0" min="0.0" help="Minimum k value."/> - <param argument="kmax" type="float" value="20" min="0.0" help="Maximum k value."/> - <param argument="kweight" type="float" value="2" help="Exponent for weighting spectra by raising k to this power."/> - <param argument="dk" type="float" value="4" help="Tapering parameter for Fourier Transform window."/> - <param argument="window" type="select" help="Fourier Transform window type."> - <option value="hanning">Hanning (cosine-squared taper)</option> - <option value="parzen">Parzen (linear taper)</option> - <option value="welch">Welch (quadratic taper)</option> - <option value="gaussian">Gaussian function window</option> - <option value="sine">Sine function window</option> - <option value="kaiser" selected="true">Kaiser-Bessel function-derived window</option> - </param> - <param argument="rmin" type="float" value="0.0" min="0.0" help="Minimum radial distance."/> - <param argument="rmax" type="float" value="10.0" min="0.0" help="Maximum radial distance."/> + <xml name="plot_limits_energy"> + <param name="x_limit_min" type="float" label="Minimum plot energy (eV)" optional="true" help="If set, plot will be limited to this value on the x axis."/> + <param name="x_limit_max" type="float" label="Maximum plot energy (eV)" optional="true" help="If set, plot will be limited to this value on the x axis."/> </xml> </macros> \ No newline at end of file