diff --git a/.codecov.yml b/.codecov.yml
new file mode 100644
index 0000000..2c5c55c
--- /dev/null
+++ b/.codecov.yml
@@ -0,0 +1,9 @@
+coverage:
+  status:
+    patch:
+      default:
+        target: 70%
+        threshold: 2%
+    project:
+      default:
+        target: auto
diff --git a/.github/workflows/clang-tidy.yml b/.github/workflows/clang-tidy.yml
index 1ef2717..f07cf06 100644
--- a/.github/workflows/clang-tidy.yml
+++ b/.github/workflows/clang-tidy.yml
@@ -24,5 +24,5 @@ jobs:
         run: |
           run-clang-tidy-19 \
             -p build \
-            -header-filter='(libdedx|tests|examples|buildbins)/' \
-            '(libdedx|tests|examples|buildbins)/.*\.c$'
+            -header-filter='(libdedx|tests|examples)/' \
+            '(libdedx|tests|examples)/.*\.c$'
diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md
index 077fb03..0faf66e 100644
--- a/CONTRIBUTING.md
+++ b/CONTRIBUTING.md
@@ -2,6 +2,30 @@
 
 ## C coding style
 
+### Symbol naming and linkage
+
+Use symbol names to make visibility obvious:
+
+- Public API functions and types use the `dedx_*` prefix.
+- Shared internal helpers used across multiple `.c` files use the `dedx_internal_*` prefix.
+- File-local helpers should be declared `static` and use short unprefixed names.
+
+Examples:
+
+```c
+/* public API */
+float dedx_get_stp(dedx_workspace *ws, dedx_config *config, float energy, int *err);
+
+/* shared internal helper */
+float dedx_internal_get_atom_mass(int id, int *err);
+
+/* file-local helper */
+static int check_energy_bounds(_dedx_lookup_data *data, float energy);
+```
+
+Do not introduce new `_dedx_*` identifiers. Existing ones should be migrated to
+the scheme above as code is touched.
+
 ### Variable declarations
 
 Declare all variables at the top of their enclosing block, before any statements.
@@ -36,10 +60,12 @@ fine when the variable is genuinely local to that scope.
 ## Thread safety
 
 libdedx is currently **not thread-safe**. There is no synchronization around
-the static path cache in `_dedx_get_data_path()`, nor around workspace mutation
-in `dedx_load_config()` / `_dedx_load_data()`. Do not share a `dedx_workspace`
+the static path cache in the file-local data-path helper in
+[`libdedx/dedx_file_access.c`](libdedx/dedx_file_access.c),
+nor around workspace mutation in `dedx_load_config()` / the internal dataset
+loading helpers. Do not share a `dedx_workspace`
 across threads without external locking.
 
-The intended fix is to replace the `done` flag in `_dedx_get_data_path()` with
-C11 `call_once()`, and audit the rest of the library for shared mutable state.
-This is tracked as a known issue.
+The intended fix is to replace the cached-path `done` flag with C11
+`call_once()`, and audit the rest of the library for shared mutable state. This
+is tracked as a known issue.
diff --git a/README.md b/README.md
index 80f9bf4..d18f8ee 100644
--- a/README.md
+++ b/README.md
@@ -1,7 +1,7 @@
 # libdedx
 
 [![CI](https://github.com/APTG/libdedx/actions/workflows/ci.yml/badge.svg)](https://github.com/APTG/libdedx/actions/workflows/ci.yml)
-[![Coverage](https://codecov.io/gh/APTG/libdedx/graph/badge.svg?branch=main)](https://codecov.io/gh/APTG/libdedx)[![License: GPL v3](https://img.shields.io/badge/License-GPLv3-blue.svg)](https://www.gnu.org/licenses/gpl-3.0)
+[![Coverage](https://codecov.io/gh/APTG/libdedx/branch/main/graph/badge.svg)](https://app.codecov.io/gh/APTG/libdedx/tree/main)[![License: GPL v3](https://img.shields.io/badge/License-GPLv3-blue.svg)](https://www.gnu.org/licenses/gpl-3.0)
 
 A C library for stopping power calculations (dE/dx) — the energy loss of a charged particle per unit length of material.
 
@@ -20,6 +20,11 @@ Full documentation: https://aptg.github.io/libdedx/
 | `DEDX_BETHE_EXT00` | Bethe formula with extensions                 | All ions      |
 | `DEDX_ICRU`     | Auto-selects ICRU49 or ICRU73 by ion type        | p, He, heavy  |
 
+Note:
+For compound targets, `libdedx` currently extends the native coverage of some original programs by falling back to Bragg/stoichiometric weighting when the upstream database does not provide that compound directly. This applies in particular to programs such as `DEDX_PSTAR`, `DEDX_ASTAR`, and `DEDX_MSTAR`, whose original target coverage is more limited than the full set of ICRU/ESTAR-style materials exposed by `libdedx`.
+
+This behavior is useful, but it means that a result returned under a given program label is not always a direct value from the original upstream program. Future releases may separate native program coverage more explicitly from these weighted extension modes.
+
 ## Quick start
 
 One-call API for a single stopping power value:
diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt
index 95cb24b..cc4f85c 100644
--- a/examples/CMakeLists.txt
+++ b/examples/CMakeLists.txt
@@ -36,4 +36,7 @@ foreach(example dedx_example dedx_bethe dedx_list dedx_use_wrappers dedx_get_sam
     add_test(NAME ${example} COMMAND ${example})
 endforeach()
 
-install (TARGETS getdedx DESTINATION bin)
\ No newline at end of file
+add_test(NAME getdedx_lookup COMMAND getdedx 2 1 276 100)
+add_test(NAME getdedx_list_programs COMMAND getdedx -1 1 276 100)
+
+install (TARGETS getdedx DESTINATION bin)
diff --git a/examples/dedx_get.c b/examples/dedx_get.c
index 1f49b57..4072fe1 100644
--- a/examples/dedx_get.c
+++ b/examples/dedx_get.c
@@ -8,6 +8,7 @@
 
 int main(int argc, char *argv[]) {
     int err = 0;
+    int exit_code = 0;
     int prog = -1;
     int target = DEDX_WATER;
     int z = -1;
@@ -16,10 +17,17 @@ int main(int argc, char *argv[]) {
     int i = 0;
     float energy = -1.0;
     float stp = 0;
-    dedx_workspace *ws;
+    dedx_workspace *ws = NULL;
     dedx_config *cfg = (dedx_config *) calloc(1, sizeof(dedx_config));
     char *str = (char *) malloc(100);
 
+    if (cfg == NULL || str == NULL) {
+        free(cfg);
+        free(str);
+        fprintf(stderr, "Allocation failure.\n");
+        return 1;
+    }
+
     if (argc != 5) {
         dedx_get_version(&vmaj, &vmin, &vpatch);
         printf("\n This is getdedx using libdEdx version %i.%i.%i", vmaj, vmin, vpatch);
@@ -35,7 +43,7 @@ int main(int argc, char *argv[]) {
         printf("   %s MSTAR CARBON WATER 100\n", argv[0]);
         printf("   %s PSTAR HYDROGEN -1 100\n", argv[0]);
         printf("\n");
-        return 0;
+        goto cleanup;
     }
 
     if (isdigit(*argv[1]))
@@ -98,7 +106,7 @@ int main(int argc, char *argv[]) {
             i++;
         }
         printf("\n");
-        return 0;
+        goto cleanup;
     }
 
     if (z == -1) {
@@ -110,7 +118,7 @@ int main(int argc, char *argv[]) {
             i++;
         }
         printf("\n");
-        return 0;
+        goto cleanup;
     }
 
     if (target == -1) {
@@ -124,7 +132,7 @@ int main(int argc, char *argv[]) {
         printf("\nAdditional materials are possible by Braggs additivity rule,\n");
         printf("as long as the constituent elements are availble in the list above.");
         printf("\n");
-        return 0;
+        goto cleanup;
     }
 
     if (energy == -1) {
@@ -133,7 +141,7 @@ int main(int argc, char *argv[]) {
                dedx_get_ion_name(z),
                dedx_get_min_energy(prog, z),
                dedx_get_max_energy(prog, z));
-        return 0;
+        goto cleanup;
     }
 
     if ((energy < dedx_get_min_energy(prog, z)) || (energy > dedx_get_max_energy(prog, z))) {
@@ -142,7 +150,8 @@ int main(int argc, char *argv[]) {
                dedx_get_ion_name(z),
                dedx_get_min_energy(prog, z),
                dedx_get_max_energy(prog, z));
-        exit(-1);
+        exit_code = 1;
+        goto cleanup;
     }
 
     if (prog == DEDX_ESTAR) {
@@ -166,7 +175,8 @@ int main(int argc, char *argv[]) {
         fprintf(stderr, "dedx_initialize, error %i:", err);
         dedx_get_error_code(str, err);
         fprintf(stderr, "  %s\n", str);
-        exit(1);
+        exit_code = 1;
+        goto cleanup;
     }
 
     /* load configuration */
@@ -179,7 +189,8 @@ int main(int argc, char *argv[]) {
         fprintf(stderr, "dedx_load_config, error %i:", err);
         dedx_get_error_code(str, err);
         fprintf(stderr, "  %s\n", str);
-        exit(1);
+        exit_code = 1;
+        goto cleanup;
     }
 
     /* get stopping power */
@@ -188,7 +199,8 @@ int main(int argc, char *argv[]) {
         fprintf(stderr, "dedx_read_energy, error %i:", err);
         dedx_get_error_code(str, err);
         fprintf(stderr, "  %s\n", str);
-        exit(1);
+        exit_code = 1;
+        goto cleanup;
     }
 
     /* check if braggs additivity rule was used */
@@ -198,9 +210,11 @@ int main(int argc, char *argv[]) {
     }
     printf("1/rho dE/dx = %6.3E MeV cm2/g\n", stp);
 
-    dedx_free_workspace(ws, &err);
+cleanup:
+    if (ws != NULL)
+        dedx_free_workspace(ws, &err);
     dedx_free_config(cfg, &err);
     free(str);
 
-    return 0;
+    return exit_code;
 }
diff --git a/examples/dedx_get_sample_csv.c b/examples/dedx_get_sample_csv.c
index f3fa185..b94935d 100644
--- a/examples/dedx_get_sample_csv.c
+++ b/examples/dedx_get_sample_csv.c
@@ -1,5 +1,6 @@
 #include <dedx.h>
 #include <stdio.h>
+#include <stdlib.h>
 
 #include "dedx_wrappers.h"
 
@@ -10,6 +11,10 @@ int main() {
     const int ions[] = {DEDX_HYDROGEN, DEDX_HELIUM, DEDX_OXYGEN, DEDX_CARBON, DEDX_URANIUM, DEDX_IRON, -1};
 
     FILE *file = fopen("result.txt", "a");
+    if (file == NULL) {
+        fprintf(stderr, "Unable to open result.txt for writing.\n");
+        return 1;
+    }
     fprintf(file, "Program, Material, Ion, Energy, Stopping Power \n");
     for (int i = 0; programs[i] != -1; i++) {
         const char *prog_name = dedx_get_program_name(programs[i]);
diff --git a/examples/dedx_use_wrappers.c b/examples/dedx_use_wrappers.c
index ff383dd..73f345f 100644
--- a/examples/dedx_use_wrappers.c
+++ b/examples/dedx_use_wrappers.c
@@ -1,5 +1,6 @@
 #include <dedx.h>
 #include <stdio.h>
+#include <stdlib.h>
 
 #include "dedx_wrappers.h"
 
diff --git a/libdedx/CMakeLists.txt b/libdedx/CMakeLists.txt
index 2a9d36c..f8aef5c 100644
--- a/libdedx/CMakeLists.txt
+++ b/libdedx/CMakeLists.txt
@@ -1,7 +1,7 @@
 # ESTAR removed from installer
-set(dedx_data data/ASTAR.bin data/astarEng.dat data/mstarEng.dat data/pstarEng.bin data/ASTAR.dat data/compos.txt data/MSTAR.dat data/gas_states.dat data/pstarEng.dat data/astarEng.bin data/effective_charge.dat data/MSTAR.bin data/mstarEng.bin data/PSTAR.bin data/PSTAR.dat data/short_names data/ICRU73.dat data/ICRU73.bin data/ICRU_PSTAR.dat data/ICRU_PSTAR.bin data/ICRU_ASTAR.dat data/ICRU_ASTAR.bin data/icru73Eng.dat data/icru73Eng.bin data/ICRU73_NEW.dat data/ICRU73_NEW.bin data/icru73_newEng.dat data/icru73_newEng.bin data/composition data/betheEng.dat data/betheEng.bin data/icru_pstarEng.dat data/icru_pstarEng.bin data/icru_astarEng.dat data/icru_astarEng.bin data/atima_compos)
+set(dedx_data data/ASTAR.bin data/astarEng.dat data/mstarEng.dat data/pstarEng.bin data/ASTAR.dat data/compos.txt data/MSTAR.dat data/gas_states.dat data/pstarEng.dat data/astarEng.bin data/effective_charge.dat data/MSTAR.bin data/mstarEng.bin data/PSTAR.bin data/PSTAR.dat data/short_names data/ICRU73.dat data/ICRU73.bin data/ICRU_PSTAR.dat data/ICRU_PSTAR.bin data/ICRU_ASTAR.dat data/ICRU_ASTAR.bin data/icru73Eng.dat data/icru73Eng.bin data/ICRU73_NEW.dat data/ICRU73_NEW.bin data/icru73_newEng.dat data/icru73_newEng.bin data/composition data/betheEng.dat data/betheEng.bin data/icru_pstarEng.dat data/icru_pstarEng.bin data/icru_astarEng.dat data/icru_astarEng.bin)
 
-set(dedx_data_bin data/ASTAR.bin data/astarEng.bin data/ICRU73.bin data/icru73Eng.bin data/ICRU73_NEW.bin data/icru73_newEng.bin data/MSTAR.bin data/mstarEng.bin data/PSTAR.bin data/ICRU_PSTAR.bin data/ICRU_ASTAR.bin data/pstarEng.bin data/betheEng.bin data/icru_astarEng.bin data/icru_pstarEng.bin data/compos.txt data/effective_charge.dat data/short_names data/composition data/gas_states.dat data/atima_compos)
+set(dedx_data_bin data/ASTAR.bin data/astarEng.bin data/ICRU73.bin data/icru73Eng.bin data/ICRU73_NEW.bin data/icru73_newEng.bin data/MSTAR.bin data/mstarEng.bin data/PSTAR.bin data/ICRU_PSTAR.bin data/ICRU_ASTAR.bin data/pstarEng.bin data/betheEng.bin data/icru_astarEng.bin data/icru_pstarEng.bin data/compos.txt data/effective_charge.dat data/short_names data/composition data/gas_states.dat)
 
 # Object library: compiled once with -fPIC so the objects feed both the
 # static and the shared library without a second compilation pass.
@@ -47,7 +47,7 @@ if(NOT WIN32)
     target_link_libraries(dedx_shared PUBLIC m)
 endif()
 
-install(FILES dedx.h dedx_error.h dedx_tools.h dedx_wrappers.h DESTINATION include)
+install(FILES dedx.h dedx_elements.h dedx_error.h dedx_tools.h dedx_wrappers.h DESTINATION include)
 
 install(TARGETS dedx dedx_shared
     RUNTIME DESTINATION bin
diff --git a/libdedx/data/atima_compos b/libdedx/data/atima_compos
deleted file mode 100644
index 4e9c573..0000000
--- a/libdedx/data/atima_compos
+++ /dev/null
@@ -1,92 +0,0 @@
-1	4.271	1.0309	1.00	19.2
-2	1.894	0.15976	1.00	41.8
-3	4.597	0.59782	1.10	40.0
-4	12.046	1.0781	1.06	63.7
-5	13.093	1.0486	1.01	76.
-6	11.364	1.00	1.03	78.
-7	3.481	1.058	1.04	82.
-8	4.302	0.93942	0.99	95.0
-9	3.522	0.74562	0.95	115.
-10	3.585	0.3424	0.90	137.
-11	2.541	0.45259	0.82	149.
-12	4.302	0.71074	0.81	156.
-13	6.023	0.90519	0.83	166.
-14	4.977	0.97411	0.88	173.
-15	3.542	0.97184	1.00	173.
-16	3.885	0.89852	0.95	180.
-17	3.22	0.70827	0.97	174.
-18	2.488	0.39816	0.99	188.
-19	1.329	0.36552	0.98	190.
-20	2.014	0.62712	0.97	191.
-21	4.015	0.81707	0.98	191.
-22	5.682	0.9943	0.97	233.
-23	7.213	1.1423	0.96	245.
-24	8.33	1.2381	0.93	257.
-25	8.15	1.1222	0.91	272.
-26	8.483	0.92705	0.9	286.
-27	8.989	1.0047	0.88	297.
-28	9.125	1.2	0.9	311.
-29	8.483	1.0661	0.9	322.
-30	6.546	0.97411	0.9	330.
-31	5.104	0.84912	0.9	334.
-32	4.428	0.95	0.85	350.
-33	4.597	1.0903	0.9	347.
-34	3.65	1.0429	0.9	348.
-35	2.562	0.49715	0.91	343.
-36	1.87	0.37755	0.92	352.
-37	1.077	0.35211	0.9	363.
-38	1.787	0.57801	0.9	366.
-39	3.041	0.77773	0.9	379.
-40	4.271	1.0207	0.9	393.
-41	5.576	1.029	0.9	417.
-42	6.407	1.2542	0.88	424.
-43	0.0	1.122	0.9	428.
-44	7.256	1.1241	0.88	441.
-45	7.256	1.0882	0.88	449.
-46	6.767	1.2709	0.9	470.
-47	5.847	1.2542	0.9	470.
-48	4.597	0.90094	0.88	469.
-49	3.836	0.74093	0.9	488.
-50	3.695	0.86054	0.9	488.
-51	3.273	0.93155	0.9	487.
-52	2.938	1.0047	0.9	485.
-53	2.343	0.55379	0.96	474.
-54	1.403	0.43289	1.2	482.
-55	.860	0.32636	0.9	488.
-56	1.544	0.5131	0.88	491.
-57	2.676	0.6950	0.88	501.
-58	2.868	0.72591	0.85	523.
-59	2.895	0.71202	0.90	535.
-60	2.923	0.67413	0.90	546.
-61	0.00	0.71418	0.92	560.
-62	3.026	0.71453	0.95	574.
-63	2.084	0.5911	0.99	580.
-64	3.026	0.70263	1.03	591.
-65	3.136	0.68049	1.05	614.
-66	3.170	0.68203	1.07	628.
-67	3.220	0.68121	1.08	650.
-68	3.273	0.68532	1.10	658.
-69	3.327	0.68715	1.08	674.
-70	2.428	0.61884	1.08	684.
-71	3.383	0.71801	1.08	694.
-72	4.428	0.83048	1.08	705.
-73	5.525	1.1222	1.09	718.
-74	6.320	1.2381	1.09	727.
-75	6.805	1.045	1.10	736.
-76	7.144	1.0733	1.11	746.
-77	7.052	1.0953	1.12	757.
-78	6.618	1.2381	1.13	790.
-79	5.904	1.2879	1.14	790.
-80	4.069	0.78654	1.15	800.
-81	3.501	0.66401	1.17	810.
-82	3.291	0.84912	1.20	823.
-83	2.827	0.88433	1.18	823.
-84	2.653	0.80746	1.17	830.
-85	0.00	0.43357	1.17	825.
-86	0.00	0.41923	1.16	794.
-87	0.00	0.43638	1.16	827.
-88	1.338	0.51464	1.16	826.
-89	0.00	0.73087	1.16	841.
-90	3.026	0.81065	1.16	847.
-91	4.015	1.9578	1.16	878.
-92	4.818	1.0257	1.16	890.
\ No newline at end of file
diff --git a/libdedx/dedx.c b/libdedx/dedx.c
index aec8a7f..310dd0c 100644
--- a/libdedx/dedx.c
+++ b/libdedx/dedx.c
@@ -17,6 +17,9 @@
 #include "dedx.h"
 
 #include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
 
 #include "dedx_bethe.h"
 #include "dedx_config.h"
@@ -29,52 +32,39 @@
 #include "dedx_spline.h"
 #include "dedx_stopping_data.h"
 #include "dedx_validate.h"
-#include "dedx_workspace.h"
 
-int _dedx_find_data(stopping_data *data, dedx_config *config, float *energy, int *err);
-
-int _dedx_check_energy_bounds(_dedx_lookup_data *data, float energy);
-int _dedx_load_bethe(stopping_data *data, int ion, int target, float *energy, int *err);
-int _dedx_load_bethe2(
-    stopping_data *data, float PZ, float PA, float TZ, float TA, float rho, float pot, float *energy, int *err);
-int _dedx_load_bethe_config(dedx_workspace *ws, float PZ, float PA, float TZ, float TA, float rho, float pot, int *err);
-int _calculate_bethe_energi_test(int ion, int target, float pot, int *err);
-int _dedx_load_data(dedx_workspace *ws, stopping_data *data, float *energy, int prog, int *err);
-
-float _dedx_get_min_energy_icru(int ion);
-float _dedx_get_max_energy_icru(int ion);
-
-int _dedx_check_ion(int prog, int ion);
-int _dedx_find_data(stopping_data *data, dedx_config *config, float *energy, int *err);
-
-int _dedx_load_config_clean(dedx_workspace *ws, dedx_config *config, int *err);
-int _dedx_load_compound(dedx_workspace *ws, dedx_config *config, int *err);
-int _dedx_load_bethe_2(stopping_data *data, dedx_config *config, float *energy, int *err);
-int _dedx_find_bragg_data_2(stopping_data *data, dedx_config *config, float *energy, int *err);
-int _dedx_load_atima(stopping_data *data, dedx_config *config, float *energy, int *err);
+static int load_data(dedx_workspace *ws, stopping_data *data, float *energy, int prog, int *err);
+static int check_energy_bounds(_dedx_lookup_data *data, float energy);
+static float get_min_energy_icru(int ion);
+static float get_max_energy_icru(int ion);
+static int check_ion(int prog, int ion);
+static int load_config_clean(dedx_workspace *ws, dedx_config *config, int *err);
+static int find_data(stopping_data *data, dedx_config *config, float *energy, int *err);
+static int load_compound(dedx_workspace *ws, dedx_config *config, int *err);
+static int load_bethe_2(stopping_data *data, dedx_config *config, float *energy, int *err);
 
 dedx_workspace *dedx_allocate_workspace(unsigned int count, int *err) {
     int i = 0;
     *err = DEDX_OK;
 
-    dedx_workspace *temp = malloc(sizeof(dedx_workspace));
+    dedx_workspace *temp = calloc(1, sizeof(dedx_workspace));
     if (temp == NULL) {
         *err = DEDX_ERR_NO_MEMORY;
         return NULL;
     }
-    temp->loaded_data = malloc(count * sizeof(_dedx_lookup_data *));
+    temp->loaded_data = (_dedx_lookup_data **) calloc(count, sizeof(_dedx_lookup_data *));
     if (temp->loaded_data == NULL) {
         *err = DEDX_ERR_NO_MEMORY;
         free(temp);
         return NULL;
     }
     for (i = 0; i < count; i++) {
-        temp->loaded_data[i] = malloc(sizeof(_dedx_lookup_data));
+        temp->loaded_data[i] = calloc(1, sizeof(_dedx_lookup_data));
         if (temp->loaded_data[i] == NULL) { /* LCOV_EXCL_START */
             int j;
             for (j = 0; j < i; j++)
                 free(temp->loaded_data[j]);
-            free(temp->loaded_data);
+            free((void *) temp->loaded_data);
             free(temp);
             *err = DEDX_ERR_NO_MEMORY;
             return NULL;
@@ -89,41 +79,16 @@ void dedx_free_workspace(dedx_workspace *workspace, int *err) {
     int i = 0;
     *err = DEDX_OK;
 
+    if (workspace == NULL)
+        return;
+
     for (i = 0; i < workspace->datasets; i++) {
         free(workspace->loaded_data[i]);
     }
-    free(workspace->loaded_data);
+    free((void *) workspace->loaded_data);
     free(workspace);
 }
 
-int _dedx_load_data(dedx_workspace *ws, stopping_data *data, float *energy, int prog, int *err) {
-    int active_dataset = ws->active_datasets;
-    *err = DEDX_OK;
-
-    _dedx_calculate_coefficients(ws->loaded_data[active_dataset]->base, energy, data->data, data->length);
-    ws->loaded_data[active_dataset]->acc.cache = 0;
-    ws->loaded_data[active_dataset]->n = data->length;
-    ws->loaded_data[active_dataset]->prog = prog;
-    ws->loaded_data[active_dataset]->ion = data->ion;
-    ws->loaded_data[active_dataset]->target = data->target;
-    ws->loaded_data[active_dataset]->datapoints = data->length;
-
-    // Increment reference Id
-    ws->active_datasets++;
-    return active_dataset;
-}
-
-/*Check the whether the energy are inside the boundary*/
-int _dedx_check_energy_bounds(_dedx_lookup_data *data, float energy) {
-    int length = data->datapoints;
-    float low = data->base[0].x;
-    float high = data->base[length - 1].x;
-    if (energy < low || energy > high) {
-        return DEDX_ERR_ENERGY_OUT_OF_RANGE;
-    }
-    return 0;
-}
-
 /*Return an explanation to the error code*/
 void dedx_get_error_code(char *err_str, int err) {
     switch (err) {
@@ -160,9 +125,6 @@ void dedx_get_error_code(char *err_str, int err) {
     case DEDX_ERR_NO_COMPOSITION:
         strcpy(err_str, "Unable to read composition file.");
         break;
-    case DEDX_ERR_NO_ATIMA_FILE:
-        strcpy(err_str, "Unable to read atima composition file.");
-        break;
     case DEDX_ERR_ENERGY_OUT_OF_RANGE:
         strcpy(err_str, "Energy out of bounds.");
         break;
@@ -230,11 +192,11 @@ void dedx_get_version(int *major, int *minor, int *patch) {
 }
 
 void dedx_get_composition(int target, float composition[][2], unsigned int *comp_len, int *err) {
-    _dedx_get_composition(target, composition, comp_len, err);
+    dedx_internal_get_composition(target, composition, comp_len, err);
 }
 
 float dedx_get_i_value(int target, int *err) {
-    return _dedx_get_i_value(target, DEDX_GAS, err);
+    return dedx_internal_get_i_value(target, DEDX_GAS, err);
 }
 
 const int *dedx_get_program_list(void) {
@@ -279,16 +241,16 @@ float dedx_get_min_energy(int program, int ion) {
         energy_min = 0.001;
         break;
     case DEDX_ICRU73_OLD:
-        energy_min = _dedx_get_min_energy_icru(ion);
+        energy_min = get_min_energy_icru(ion);
         break;
     case DEDX_ICRU73:
-        energy_min = _dedx_get_min_energy_icru(ion);
+        energy_min = get_min_energy_icru(ion);
         break;
     case DEDX_ICRU49:
-        energy_min = _dedx_get_min_energy_icru(ion);
+        energy_min = get_min_energy_icru(ion);
         break;
     case DEDX_ICRU:
-        energy_min = _dedx_get_min_energy_icru(ion);
+        energy_min = get_min_energy_icru(ion);
         break;
     case DEDX_DEFAULT:
         energy_min = 0.001;
@@ -317,16 +279,16 @@ float dedx_get_max_energy(int program, int ion) {
         energy_max = 1000.0;
         break;
     case DEDX_ICRU73_OLD:
-        energy_max = _dedx_get_max_energy_icru(ion);
+        energy_max = get_max_energy_icru(ion);
         break;
     case DEDX_ICRU73:
-        energy_max = _dedx_get_max_energy_icru(ion);
+        energy_max = get_max_energy_icru(ion);
         break;
     case DEDX_ICRU49:
-        energy_max = _dedx_get_max_energy_icru(ion);
+        energy_max = get_max_energy_icru(ion);
         break;
     case DEDX_ICRU:
-        energy_max = _dedx_get_max_energy_icru(ion);
+        energy_max = get_max_energy_icru(ion);
         break;
     case DEDX_DEFAULT:
         energy_max = 1000.0;
@@ -338,48 +300,192 @@ float dedx_get_max_energy(int program, int ion) {
     return energy_max;
 }
 
-float _dedx_get_min_energy_icru(int ion) {
-    float _energy_min = 0;
+int dedx_load_config(dedx_workspace *ws, dedx_config *config, int *err) {
+    int cfg_id = -1;
+
+    config->loaded = 0;
+    config->cfg_id = -1;
+    dedx_internal_validate_config(config, err);
+    if (*err != 0)
+        return -1;
+    if (config->elements_id != NULL)
+        cfg_id = load_compound(ws, config, err);
+    else
+        cfg_id = load_config_clean(ws, config, err);
+    if (*err != 0 || cfg_id < 0)
+        return -1;
+    dedx_internal_set_names(config, err);
+    if (*err != 0)
+        return -1;
+    config->cfg_id = cfg_id;
+    config->loaded = 1;
+    return 0;
+}
+
+float dedx_get_stp(dedx_workspace *ws, dedx_config *config, float energy, int *err) {
+    int id = config->cfg_id;
+
+    *err = DEDX_OK;
+
+    if (id < 0 || id >= ws->active_datasets) {
+        *err = DEDX_ERR_INVALID_DATASET_ID;
+        return 0;
+    }
+
+    // Check that the energy is within the boundary
+    *err = check_energy_bounds(ws->loaded_data[id], energy);
+    if (*err != DEDX_OK) {
+        return 0;
+    }
+
+    // Evaluating the spline function
+    return dedx_internal_evaluate_spline(
+        ws->loaded_data[id]->base, energy, &(ws->loaded_data[id]->acc), ws->loaded_data[id]->n);
+}
+
+void dedx_free_config(dedx_config *config, int *err) {
+    if (config != NULL) {
+        if (config->elements_id != NULL)
+            free(config->elements_id);
+        if (config->elements_atoms != NULL)
+            free(config->elements_atoms);
+        if (config->elements_mass_fraction != NULL)
+            free(config->elements_mass_fraction);
+        if (config->elements_i_value != NULL)
+            free(config->elements_i_value);
+        free(config);
+    }
+    *err = DEDX_OK;
+}
+
+float dedx_get_simple_stp(int ion, int target, float energy, int *err) {
+    dedx_config *config = NULL;
+    dedx_workspace *ws = NULL;
+    int cleanup_err = DEDX_OK;
+    float stp;
+
+    *err = DEDX_OK;
+
+    config = calloc(1, sizeof(dedx_config));
+    if (config == NULL) {
+        *err = DEDX_ERR_NO_MEMORY;
+        return 0.0f;
+    }
+
+    config->ion = ion;
+    config->target = target;
+    config->program = DEDX_ICRU;
+
+    ws = dedx_allocate_workspace(1, err);
+    if (*err != DEDX_OK) {
+        dedx_free_config(config, &cleanup_err);
+        return 0.0f;
+    }
+
+    if (dedx_load_config(ws, config, err) != 0) {
+        dedx_free_config(config, &cleanup_err);
+        config = calloc(1, sizeof(dedx_config));
+        if (config == NULL) {
+            *err = DEDX_ERR_NO_MEMORY;
+            dedx_free_workspace(ws, &cleanup_err);
+            return 0.0f;
+        }
+
+        config->ion = ion;
+        config->target = target;
+        config->program = DEDX_DEFAULT;
+        if (dedx_load_config(ws, config, err) != 0) {
+            dedx_free_config(config, &cleanup_err);
+            dedx_free_workspace(ws, &cleanup_err);
+            return 0.0f;
+        }
+    }
+
+    stp = dedx_get_stp(ws, config, energy, err);
+    dedx_free_config(config, &cleanup_err);
+    dedx_free_workspace(ws, &cleanup_err);
+    if (*err != DEDX_OK)
+        return 0.0f;
+
+    return stp;
+}
+
+static int load_data(dedx_workspace *ws, stopping_data *data, float *energy, int prog, int *err) {
+    int active_dataset = ws->active_datasets;
+
+    *err = DEDX_OK;
+
+    if (active_dataset >= ws->datasets) {
+        *err = DEDX_ERR_INVALID_DATASET_ID;
+        return -1;
+    }
+
+    dedx_internal_calculate_coefficients(ws->loaded_data[active_dataset]->base, energy, data->data, data->length);
+
+    ws->loaded_data[active_dataset]->acc.cache = 0;
+    ws->loaded_data[active_dataset]->n = data->length;
+    ws->loaded_data[active_dataset]->prog = prog;
+    ws->loaded_data[active_dataset]->ion = data->ion;
+    ws->loaded_data[active_dataset]->target = data->target;
+    ws->loaded_data[active_dataset]->datapoints = data->length;
+
+    ws->active_datasets++;
+    return active_dataset;
+}
+
+static int check_energy_bounds(_dedx_lookup_data *data, float energy) {
+    int length = data->datapoints;
+    float low = data->base[0].x;
+    float high = data->base[length - 1].x;
+
+    if (energy < low || energy > high) {
+        return DEDX_ERR_ENERGY_OUT_OF_RANGE;
+    }
+    return 0;
+}
+
+static float get_min_energy_icru(int ion) {
+    float energy_min = 0;
+
     switch (ion) {
     case 1:
-        _energy_min = 0.001;
+        energy_min = 0.001;
         break;
     case 2:
-        _energy_min = 0.001 / 4.0;
+        energy_min = 0.001 / 4.0;
         break;
     default:
-        _energy_min = 0.025;
+        energy_min = 0.025;
         break;
     }
-    return _energy_min;
+    return energy_min;
 }
 
-float _dedx_get_max_energy_icru(int ion) {
-    float _energy_max = 0;
+static float get_max_energy_icru(int ion) {
+    float energy_max = 0;
+
     switch (ion) {
     case 1:
-        _energy_max = 10000.0;
+        energy_max = 10000.0;
         break;
     case 2:
-        _energy_max = 1000.0 / 4.0;
+        energy_max = 1000.0 / 4.0;
         break;
     default:
-        _energy_max = 1000.0;
+        energy_max = 1000.0;
         break;
     }
-    return _energy_max;
+    return energy_max;
 }
 
-int _dedx_check_ion(int prog, int ion) {
-    // checks if ion is available in program. Returns 0 if yes, else -1.
+static int check_ion(int prog, int ion) {
     const int *ion_list;
     int i = 0;
 
     if (prog >= DEDX_DEFAULT) {
         if ((ion < 1) || (ion > 120))
             return 0;
-        else
-            return 1;
+        return 1;
     }
 
     ion_list = dedx_get_ion_list(prog);
@@ -391,45 +497,28 @@ int _dedx_check_ion(int prog, int ion) {
     return 0;
 }
 
-void dedx_load_config(dedx_workspace *ws, dedx_config *config, int *err) {
-    int cfg_id;
-    _dedx_validate_config(config, err);
-    if (*err != 0)
-        return;
-    if (config->elements_id != NULL)
-        cfg_id = _dedx_load_compound(ws, config, err);
-    else
-        cfg_id = _dedx_load_config_clean(ws, config, err);
-    _dedx_set_names(config, err);
-    config->cfg_id = cfg_id;
-    config->loaded = 1;
-}
-
-int _dedx_load_config_clean(dedx_workspace *ws, dedx_config *config, int *err) {
-    float energy[_DEDX_MAXELEMENTS];
+static int load_config_clean(dedx_workspace *ws, dedx_config *config, int *err) {
+    float energy[DEDX_MAX_ELEMENTS];
     int cfg;
     int prog = config->program;
     int ion = config->ion;
     int target = config->target;
+    stopping_data data;
+
     config->bragg_used = 0;
     *err = DEDX_OK;
-    stopping_data data;
 
-    // check if ion is available in requested program
-    if (!_dedx_check_ion(prog, ion)) {
+    if (!check_ion(prog, ion)) {
         *err = DEDX_ERR_ION_NOT_SUPPORTED;
         return -1;
     }
     config->_temp_i_value = config->i_value;
-    // Load data
-    _dedx_find_data(&data, config, energy, err);
+    find_data(&data, config, energy, err);
 
     if (*err != 0) {
-        // Check whether the error was that the combination wasn't in the data files and target is a compound
-
         if (*err == DEDX_ERR_COMBINATION_NOT_FOUND && target > 99) {
             *err = DEDX_OK;
-            _dedx_evaluate_compound(config, err);
+            dedx_internal_evaluate_compound(config, err);
             if (*err != 0)
                 return -1;
             if (config->elements_length == 0) {
@@ -437,23 +526,22 @@ int _dedx_load_config_clean(dedx_workspace *ws, dedx_config *config, int *err) {
                 *err = DEDX_ERR_TARGET_NOT_FOUND;
                 return -1;
             }
-            cfg = _dedx_load_compound(ws, config, err);
+            cfg = load_compound(ws, config, err);
             if (*err != 0)
                 return -1;
         }
-        if (*err != 0) {
+        if (*err != 0)
             return -1;
-        }
-    } else
-        cfg = _dedx_load_data(ws, &data, energy, prog, err);
+    } else {
+        cfg = load_data(ws, &data, energy, prog, err);
+    }
     return cfg;
 }
 
-int _dedx_find_data(stopping_data *data, dedx_config *config, float *energy, int *err) {
+static int find_data(stopping_data *data, dedx_config *config, float *energy, int *err) {
     int prog = config->program;
     int target = config->target;
     int ion = config->ion;
-
     int prog_load = prog;
     int ion_load = ion;
     int target_load = target;
@@ -465,12 +553,10 @@ int _dedx_find_data(stopping_data *data, dedx_config *config, float *energy, int
             prog_load = _DEDX_0008;
         } else if (ion == 2) {
             prog_load = DEDX_ICRU49;
+        } else if (!(target == DEDX_WATER || target == DEDX_AIR)) {
+            prog_load = DEDX_ICRU73_OLD;
         } else {
-            if (!(target == DEDX_WATER || target == DEDX_AIR)) {
-                prog_load = DEDX_ICRU73_OLD;
-            } else {
-                prog_load = DEDX_ICRU73;
-            }
+            prog_load = DEDX_ICRU73;
         }
     } else if (prog == DEDX_ICRU49) {
         if (ion == 1) {
@@ -481,52 +567,53 @@ int _dedx_find_data(stopping_data *data, dedx_config *config, float *energy, int
             *err = DEDX_ERR_COMBINATION_NOT_FOUND;
             return -1;
         }
-    }
-
-    // ESTAR not supported
-    else if (prog == DEDX_ESTAR) {
+    } else if (prog == DEDX_ESTAR) {
         *err = DEDX_ERR_ESTAR_NOT_IMPL;
         return -1;
     } else if (prog == DEDX_ICRU73 && !(target == 276 || target == 104)) {
         prog_load = DEDX_ICRU73_OLD;
-    } else if (prog == DEDX_MSTAR && ion > 1) // Load ASTAR data, and scale them with the MSTAR method
-    {
+    } else if (prog == DEDX_MSTAR && ion > 1) {
         ion_load = 2;
-    } else if (prog == DEDX_BETHE_EXT00 || prog == DEDX_DEFAULT) // Load bethe data
-    {
-        prog_load = 101;
+    } else if (prog == DEDX_BETHE_EXT00 || prog == DEDX_DEFAULT) {
         int prog_temp = config->program;
-        _dedx_read_energy_data(energy, prog_load, err);
+
+        prog_load = 101;
+        dedx_internal_read_energy_data(energy, prog_load, err);
         config->program = prog_load;
-        _dedx_load_bethe_2(data, config, energy, err);
+        load_bethe_2(data, config, energy, err);
         config->program = prog_temp;
         return 0;
     }
-    _dedx_read_binary_data(data, prog_load, ion_load, target_load, err);
+
+    dedx_internal_read_binary_data(data, prog_load, ion_load, target_load, err);
     if (*err != 0)
         return -1;
-    _dedx_read_energy_data(energy, prog_load, err);
+
+    dedx_internal_read_energy_data(energy, prog_load, err);
     if (prog == DEDX_MSTAR) {
         stopping_data out;
         char mode = 'b';
+
         if (config->mstar_mode != '\0') {
             mode = config->mstar_mode;
         }
         data->ion = config->ion;
-        _dedx_convert_energy_to_mstar(data, &out, mode, config, energy);
+        dedx_internal_convert_energy_to_mstar(data, &out, mode, config, energy, err);
+        if (*err != 0)
+            return -1;
         memcpy(data, &out, sizeof(stopping_data));
         data->ion = ion;
     }
     return 0;
 }
 
-int _dedx_load_compound(dedx_workspace *ws, dedx_config *config, int *err) {
+static int load_compound(dedx_workspace *ws, dedx_config *config, int *err) {
     int i = 0;
     int j = 0;
     int length = config->elements_length;
     int *targets = config->elements_id;
     float *weight;
-    float energy[_DEDX_MAXELEMENTS];
+    float energy[DEDX_MAX_ELEMENTS];
     float i_value;
     int target;
     stopping_data data;
@@ -539,7 +626,6 @@ int _dedx_load_compound(dedx_workspace *ws, dedx_config *config, int *err) {
         return -1;
     }
     weight = config->elements_mass_fraction;
-    // Backup i_value and target
     i_value = config->i_value;
     target = config->target;
     for (i = 0; i < length; i++) {
@@ -548,16 +634,17 @@ int _dedx_load_compound(dedx_workspace *ws, dedx_config *config, int *err) {
             config->_temp_i_value = config->elements_i_value[i];
             if (config->elements_i_value[i] <= 0.0) {
                 *err = DEDX_ERR_INVALID_I_VALUE;
+                free(compound_data);
                 return -1;
             }
         }
-        _dedx_find_data(&compound_data[i], config, energy, err);
+        find_data(&compound_data[i], config, energy, err);
         if (*err != 0) {
             free(compound_data);
             return -1;
         }
     }
-    // Restore I_value and target
+
     config->i_value = i_value;
     config->target = target;
     for (j = 0; j < compound_data[0].length; j++) {
@@ -568,32 +655,33 @@ int _dedx_load_compound(dedx_workspace *ws, dedx_config *config, int *err) {
     }
     data.length = compound_data[0].length;
     free(compound_data);
-    return _dedx_load_data(ws, &data, energy, config->program, err);
+    return load_data(ws, &data, energy, config->program, err);
 }
 
-int _dedx_load_bethe_2(stopping_data *data, dedx_config *config, float *energy, int *err) {
+static int load_bethe_2(stopping_data *data, dedx_config *config, float *energy, int *err) {
     int i = 0;
+    float PZ, PA, TZ, TA, rho, pot;
+
     *err = DEDX_OK;
     if (config->target > 99) {
         *err = DEDX_ERR_COMBINATION_NOT_FOUND;
         return -1;
     }
-    float PZ, PA, TZ, TA, rho, pot;
-    PZ = _dedx_get_atom_charge(config->ion, err);
-    PA = _dedx_get_atom_mass(config->ion, err);
-    TZ = _dedx_get_atom_charge(config->target, err);
-    TA = _dedx_get_atom_mass(config->target, err);
+
+    PZ = dedx_internal_get_atom_charge(config->ion, err);
+    PA = dedx_internal_get_atom_mass(config->ion, err);
+    TZ = dedx_internal_get_atom_charge(config->target, err);
+    TA = dedx_internal_get_atom_mass(config->target, err);
     rho = config->rho;
     pot = config->_temp_i_value;
     data->length = 122;
-    // Get energy grid.
-    _dedx_read_energy_data(energy, DEDX_BETHE_EXT00, err);
+    dedx_internal_read_energy_data(energy, DEDX_BETHE_EXT00, err);
     if (*err != 0)
-        return 0;
-    // Fill the data grid with values.
+        return -1;
+
     _dedx_bethe_coll_struct *bethe = (_dedx_bethe_coll_struct *) calloc(1, sizeof(_dedx_bethe_coll_struct));
     for (i = 0; i < data->length; i++) {
-        data->data[i] = _dedx_calculate_bethe_energy(bethe, energy[i], PZ, PA, TZ, TA, rho, pot);
+        data->data[i] = dedx_internal_calculate_bethe_energy(bethe, energy[i], PZ, PA, TZ, TA, rho, pot);
     }
     if (bethe->bet != NULL)
         free(bethe->bet);
@@ -602,67 +690,3 @@ int _dedx_load_bethe_2(stopping_data *data, dedx_config *config, float *energy,
     free(bethe);
     return 0;
 }
-
-int _dedx_load_atima(stopping_data *data, dedx_config *config, float *energy, int *err) {
-    *err = DEDX_OK;
-    return 0;
-}
-
-float dedx_get_stp(dedx_workspace *ws, dedx_config *config, float energy, int *err) {
-    int id = config->cfg_id;
-    *err = DEDX_OK;
-    // Check that the energy is within the boundary
-    if ((*err = _dedx_check_energy_bounds(ws->loaded_data[id], energy)) != DEDX_OK)
-        return 0;
-    if (id > ws->active_datasets) // Check that the dataset is loaded.
-    {
-        *err = DEDX_ERR_INVALID_DATASET_ID;
-        return 0;
-    }
-    // Evaluating the spline function
-    return _dedx_evaluate_spline(
-        ws->loaded_data[id]->base, energy, &(ws->loaded_data[id]->acc), ws->loaded_data[id]->n);
-}
-
-void dedx_free_config(dedx_config *config, int *err) {
-    if (config != NULL) {
-        if (config->elements_id != NULL)
-            free(config->elements_id);
-        if (config->elements_atoms != NULL)
-            free(config->elements_atoms);
-        if (config->elements_mass_fraction != NULL)
-            free(config->elements_mass_fraction);
-        if (config->elements_i_value != NULL)
-            free(config->elements_i_value);
-        free(config);
-    }
-    *err = DEDX_OK;
-}
-
-float dedx_get_simple_stp(int ion, int target, float energy, int *err) {
-    dedx_config *config = (dedx_config *) calloc(1, sizeof(dedx_config));
-    float stp;
-    config->target = target;
-    config->ion = ion;
-    config->program = DEDX_ICRU;
-    dedx_workspace *ws = dedx_allocate_workspace(1, err);
-    if (*err != 0)
-        return 0;
-    dedx_load_config(ws, config, err);
-    if (*err != 0) {
-        dedx_free_config(config, err);
-        config = (dedx_config *) calloc(1, sizeof(dedx_config));
-        config->ion = ion;
-        config->target = target;
-        config->program = 100;
-        dedx_load_config(ws, config, err);
-        if (*err != 0)
-            return 0;
-    }
-    stp = dedx_get_stp(ws, config, energy, err);
-    if (*err != 0)
-        return 0;
-    dedx_free_config(config, err);
-    dedx_free_workspace(ws, err);
-    return stp;
-}
diff --git a/libdedx/dedx.h b/libdedx/dedx.h
index 14ec68b..7eed689 100644
--- a/libdedx/dedx.h
+++ b/libdedx/dedx.h
@@ -1,6 +1,7 @@
-#ifndef DEDX_H_INCLUDED
-#define DEDX_H_INCLUDED
+#ifndef DEDX_H
+#define DEDX_H
 
+#include "dedx_elements.h"
 #include "dedx_error.h"
 
 /**
@@ -25,8 +26,6 @@
  * @endcode
  */
 
-#define _DEDX_MAXELEMENTS 150
-
 /**
  * @defgroup programs Stopping power programs
  * @brief Identifiers for the supported stopping power databases/models.
@@ -67,339 +66,17 @@ enum {
  * @{
  */
 enum {
-    DEDX_MSTAR_MODE_A = 'a',
-    DEDX_MSTAR_MODE_B = 'b',
-    DEDX_MSTAR_MODE_G = 'g',
-    DEDX_MSTAR_MODE_H = 'h',
-    DEDX_MSTAR_MODE_C = 'c',
-    DEDX_MSTAR_MODE_D = 'd',
+    DEDX_MSTAR_MODE_A = 'a',      /**< Automatic state selection; use condensed/gas base mode (`c` or `g`). */
+    DEDX_MSTAR_MODE_B = 'b',      /**< Automatic state selection; prefer special state mode (`d` or `h`). */
+    DEDX_MSTAR_MODE_G = 'g',      /**< Gaseous target mode. */
+    DEDX_MSTAR_MODE_H = 'h',      /**< Special gaseous-target mode for supported ions. */
+    DEDX_MSTAR_MODE_C = 'c',      /**< Condensed target mode. */
+    DEDX_MSTAR_MODE_D = 'd',      /**< Special condensed-target mode. */
     DEDX_MSTAR_MODE_DEFAULT = 'b' /**< Recommended by Helmut Paul */
 };
 
 /** @} */
 
-/**
- * @defgroup ions_and_materials Ion and material identifiers
- * @brief Identifiers for projectile ions and target materials.
- *
- * Elemental ions (Z=1–98) correspond to their atomic number.
- * Compound materials follow sequentially after the elements.
- * Both share a single enum to preserve the numeric ID mapping.
- * @{
- */
-enum {
-    /* --- Elemental ions (Z = atomic number) --- */
-    DEDX_HYDROGEN = 1,
-    DEDX_HELIUM,
-    DEDX_LITHIUM,
-    DEDX_BERYLLIUM,
-    DEDX_BORON,
-    DEDX_CARBON,
-    DEDX_GRAPHITE = 906,
-    DEDX_NITROGEN = 7,
-    DEDX_OXYGEN,
-    DEDX_FLUORINE,
-    DEDX_NEON,
-    DEDX_SODIUM,
-    DEDX_MAGNESIUM,
-    DEDX_ALUMINUM,
-    DEDX_SILICON,
-    DEDX_PHOSPHORUS,
-    DEDX_SULFUR,
-    DEDX_CHLORINE,
-    DEDX_ARGON,
-    DEDX_POTASSIUM,
-    DEDX_CALCIUM,
-    DEDX_SCANDIUM,
-    DEDX_TITeANIUM,
-    DEDX_VANADIUM,
-    DEDX_CHROMIUM,
-    DEDX_MANGANESE,
-    DEDX_IRON,
-    DEDX_COBALT,
-    DEDX_NICKEL,
-    DEDX_COPPER,
-    DEDX_ZINC,
-    DEDX_GALLIUM,
-    DEDX_GERMANIUM,
-    DEDX_ARSENIC,
-    DEDX_SELENIUM,
-    DEDX_BROMINE,
-    DEDX_KRYPTON,
-    DEDX_RUBIDIUM,
-    DEDX_STRONTIUM,
-    DEDX_YTTRIUM,
-    DEDX_ZIRCONIUM,
-    DEDX_NIOBIUM,
-    DEDX_MOLYBDENUM,
-    DEDX_TECHNETIUM,
-    DEDX_RUTHENIUM,
-    DEDX_RHODIUM,
-    DEDX_PALLADIUM,
-    DEDX_SILVER,
-    DEDX_CADMIUM,
-    DEDX_INDIUM,
-    DEDX_TIN,
-    DEDX_ANTIMONY,
-    DEDX_TELLURIUM,
-    DEDX_IODINE,
-    DEDX_XENON,
-    DEDX_CESIUM,
-    DEDX_BARIUM,
-    DEDX_LANTHANUM,
-    DEDX_CERIUM,
-    DEDX_PRASEODYMIUM,
-    DEDX_NEODYMIUM,
-    DEDX_PROMETHIUM,
-    DEDX_SAMARIUM,
-    DEDX_EUROPIUM,
-    DEDX_GADOLINIUM,
-    DEDX_TERBIUM,
-    DEDX_DYSPROSIUM,
-    DEDX_HOLMIUM,
-    DEDX_ERBIUM,
-    DEDX_THULIUM,
-    DEDX_YTTERBIUM,
-    DEDX_LUTETIUM,
-    DEDX_HAFNIUM,
-    DEDX_TANTALUM,
-    DEDX_TUNGSTEN,
-    DEDX_RHENIUM,
-    DEDX_OSMIUM,
-    DEDX_IRIDIUM,
-    DEDX_PLATINUM,
-    DEDX_GOLD,
-    DEDX_MERCURY,
-    DEDX_THALLIUM,
-    DEDX_LEAD,
-    DEDX_BISMUTH,
-    DEDX_POLONIUM,
-    DEDX_ASTATINE,
-    DEDX_RADON,
-    DEDX_FRANCIUM,
-    DEDX_RADIUM,
-    DEDX_ACTINIUM,
-    DEDX_THORIUM,
-    DEDX_PROTACTINIUM,
-    DEDX_URANIUM,
-    DEDX_NEPTUNIUM,
-    DEDX_PLUTONIUM,
-    DEDX_AMERICIUM,
-    DEDX_CURIUM,
-    DEDX_BERKELIUM,
-    DEDX_CALIFORNIUM,
-
-    /* --- Compound and mixture targets (NIST material database) --- */
-    DEDX_A150_TISSUE_EQUIVALENT_PLASTIC,
-    DEDX_ACETONE,
-    DEDX_ACETYLENE,
-    DEDX_ADENINE,
-    DEDX_ADIPOSE_TISSUE_ICRP,
-    DEDX_AIR_DRY_NEAR_SEA_LEVEL,
-    DEDX_ALANINE,
-    DEDX_ALUMINUMOXIDE,
-    DEDX_AMBER,
-    DEDX_AMMONIA,
-    DEDX_ANILINE,
-    DEDX_ANTHRACENE,
-    DEDX_B100,
-    DEDX_BAKELITE,
-    DEDX_BARIUM_FLUORIDE,
-    DEDX_BARIUM_SULFATE,
-    DEDX_BENZENE,
-    DEDX_BERYLLIUM_OXIDE,
-    DEDX_BISMUTH_GERMANIUM_OXIDE,
-    DEDX_BLOOD_ICRP,
-    DEDX_BONE_COMPACT_ICRU,
-    DEDX_BONE_CORTICAL_ICRP,
-    DEDX_BORON_CARBIDE,
-    DEDX_BORON_OXIDE,
-    DEDX_BRAIN_ICRP,
-    DEDX_BUTANE,
-    DEDX_N_BUTYLALCOHOL,
-    DEDX_C552_AIR_EQUIVALENT_PLASTIC,
-    DEDX_CADMIUM_TELLURIDE,
-    DEDX_CADMIUM_TUNGSTATE,
-    DEDX_CALCIUM_CARBONATE,
-    DEDX_CALCIUM_FLUORIDE,
-    DEDX_CALCIUM_OXIDE,
-    DEDX_CALCIUM_SULFATE,
-    DEDX_CALCIUM_TUNGSTATE,
-    DEDX_CARBON_DIOXIDE,
-    DEDX_CARBON_TETRACHLORIDE,
-    DEDX_CELLULOSE_ACETATE_CELLOPHANE,
-    DEDX_CELLULOSE_ACETATE_BUTYRATE,
-    DEDX_CELLULOSE_NITRATE,
-    DEDX_CERIC_SULFATE_DOSIMETER_SOLUTION,
-    DEDX_CESIUM_FLUORIDE,
-    DEDX_CESIUM_IODIDE,
-    DEDX_CHLORO_BENZENE,
-    DEDX_CHLOROFORM,
-    DEDX_CONCRETE_PORTLAND,
-    DEDX_CYCLOHEXANE,
-    DEDX_DICHLOROBENZENE,
-    DEDX_DICHLORODIETHYL_ETHER,
-    DEDX_DICHLOROETHANE,
-    DEDX_DIETHYLETHER,
-    DEDX_N_N_DIMETHYL_FORMAMIDE,
-    DEDX_DIMETHYL_SULFOXIDE,
-    DEDX_ETHANE,
-    DEDX_ETHYL_ALCOHOL,
-    DEDX_ETHYL_CELLULOSE,
-    DEDX_ETHYLENE,
-    DEDX_EYE_LENS_ICRP,
-    DEDX_FERRIC_OXIDE,
-    DEDX_FERRO_BORIDE,
-    DEDX_FERROUS_OXIDE,
-    DEDX_FERROUS_SULFATE_DOSIMETER_SOLUTION,
-    DEDX_FREON_12,
-    DEDX_FREON_12B2,
-    DEDX_FREON_13,
-    DEDX_FREON_13B1,
-    DEDX_FREON_13I1,
-    DEDX_GADOLINIUM_OXYSULFIDE,
-    DEDX_GALLIUM_ARSENIDE,
-    DEDX_GEL_IN_PHOTOGRAPHIC_EMULSION,
-    DEDX_GLASS_PYREX,
-    DEDX_GLASS_LEAD,
-    DEDX_GLASS_PLATE, /* 169,170,171 */
-    DEDX_GLUCOSE,
-    DEDX_GLUTAMINE,
-    DEDX_GLYCEROL,
-    DEDX_GUANINE,
-    DEDX_GYPSUM_PLASTER_OF_PARIS,
-    DEDX_N_HEPTANE,
-    DEDX_N_HEXANE,
-    DEDX_KAPTON_POLYIMIDE_FILM,
-    DEDX_LANTHANUM_OXYBROMIDE,
-    DEDX_LANTHANUM_OXYSULFIDE,
-    DEDX_LEAD_OXIDE,
-    DEDX_LITHIUM_AMIDE,
-    DEDX_LITHIUM_CARBONATE,
-    DEDX_LITHIUM_FLUORIDE,
-    DEDX_LITHIUM_HYDRIDE,
-    DEDX_LITHIUM_IODIDE,
-    DEDX_LITHIUM_OXIDE,
-    DEDX_LITHIUM_TETRABORATE,
-    DEDX_LUNG_ICRP,
-    DEDX_M3_WAX,
-    DEDX_MAGNESIUM_CARBONATE,
-    DEDX_MAGNESIUM_FLUORIDE,
-    DEDX_MAGNESIUM_OXIDE,
-    DEDX_MAGNESIUM_TETRABORATE,
-    DEDX_MERCURIC_IODIDE,
-    DEDX_METHANE,
-    DEDX_METHANOL,
-    DEDX_MIX_D_WAX,
-    DEDX_MS20_TISSUE_SUBSTITUTE,
-    DEDX_MUSCLE_SKELETAL,
-    DEDX_MUSCLE_STRIATED,
-    DEDX_MUSCLE_EQUIVALENT_LIQUID_WITH_SUCROSE,
-    DEDX_MUSCLE_EQUIVALENT_LIQUID_WITHOUT_SUCROSE,
-    DEDX_NAPHTHALENE,
-    DEDX_NITROBENZENE,
-    DEDX_NITROUS_OXIDE,
-    DEDX_NYLON_DUPONT_ELVAMIDE_8062,
-    DEDX_NYLON_TYPE_6_AND_6_6,
-    DEDX_NYLON_TYPE_6_10,
-    DEDX_NYLON_TYPE_11_RILSAN,
-    DEDX_OCTANE_LIQUID,
-    DEDX_PARAFFIN_WAX,
-    DEDX_N_PENTANE,
-    DEDX_PHOTOGRAPHIC_EMULSION,
-    DEDX_PLASTIC_SCINTILLATOR_VINYLTOLUENE_BASED,
-    DEDX_PLUTONIUM_DIOXIDE,
-    DEDX_POLYACRYLONITRILE,
-    DEDX_POLYCARBONATE_MAKROLON_LEXAN,
-    DEDX_POLYCHLOROSTYRENE,
-    DEDX_POLYETHYLENE,
-    DEDX_MYLAR,
-    DEDX_LUCITE_PERSPEX_PMMA,
-    DEDX_POLYOXYMETHYLENE,
-    DEDX_POLYPROPYLENE,
-    DEDX_POLYSTYRENE,
-    DEDX_POLYTETRAFLUOROETHYLENE,
-    DEDX_POLYTRIFLUOROCHLOROETHYLENE,
-    DEDX_POLYVINYL_ACETATE,
-    DEDX_POLYVINYL_ALCOHOL,
-    DEDX_POLYVINYL_BUTYRAL,
-    DEDX_POLYVINYL_CHLORIDE,
-    DEDX_POLYVINYLIDENE_CHLORIDE_SARAN,
-    DEDX_POLYVINYLIDENE_FLUORIDE,
-    DEDX_POLYVINYL_PYRROLIDONE,
-    DEDX_POTASSIUM_IODIDE,
-    DEDX_POTASSIUM_OXIDE,
-    DEDX_PROPANE,
-    DEDX_PROPANE_LIQUID,
-    DEDX_N_PROPYL_ALCOHOL,
-    DEDX_PYRIDINE,
-    DEDX_RUBBER_BUTYL,
-    DEDX_RUBBER_NATURAL,
-    DEDX_RUBBER_NEOPRENE,
-    DEDX_SILICON_DIOXIDE,
-    DEDX_SILVER_BROMIDE,
-    DEDX_SILVER_CHLORIDE,
-    DEDX_SILVER_HALIDES_IN_PHOTOGRAPHIC_EMULSION,
-    DEDX_SILVER_IODIDE,
-    DEDX_SKIN_ICRP,
-    DEDX_SODIUM_CARBONATE,
-    DEDX_SODIUM_IODIDE,
-    DEDX_SODIUM_MONOXIDE,
-    DEDX_SODIUM_NITRATE,
-    DEDX_STILBENE,
-    DEDX_SUCROSE,
-    DEDX_TERPHENYL,
-    DEDX_TESTES_ICRP,
-    DEDX_TETRACHLOROETHYLENE,
-    DEDX_THALLIUM_CHLORIDE,
-    DEDX_TISSUE_SOFT_ICRP,
-    DEDX_TISSUE_SOFT_ICRU_FOUR_COMPONENT,
-    DEDX_TISSUE_EQUIVALENT_GAS_METHANE_BASED,
-    DEDX_TISSUE_EQUIVALENT_GAS_PROPANE_BASED,
-    DEDX_TITANIUM_DIOXIDE,
-    DEDX_TOLUENE,
-    DEDX_TRICHLOROETHYLENE,
-    DEDX_TRIETHYL_PHOSPHATE,
-    DEDX_TUNGSTEN_HEXAFLUORIDE,
-    DEDX_URANIUM_DICARBIDE,
-    DEDX_URANIUM_MONOCARBIDE,
-    DEDX_URANIUM_OXIDE,
-    DEDX_UREA,
-    DEDX_VALINE,
-    DEDX_VITON_FLUOROELASTOMER,
-    DEDX_WATER_LIQUID,
-    DEDX_WATER_VAPOR,
-    DEDX_XYLENE
-};
-
-/** @} */ /* ions_and_materials */
-
-/** @defgroup special_ions Special particle identifiers
- *  @brief Identifiers for particles that do not map directly to atomic number.
- *  @{
- */
-#define DEDX_PROTON 1        /**< Proton (Z=1, same as DEDX_HYDROGEN) */
-#define DEDX_ELECTRON 1001   /**< Electron */
-#define DEDX_POSITRON 1002   /**< Positron */
-#define DEDX_PIMINUS 1003    /**< Pi minus */
-#define DEDX_PIPLUS 1004     /**< Pi plus */
-#define DEDX_PIZERO 1005     /**< Pi zero */
-#define DEDX_ANTIPROTON 1006 /**< Antiproton */
-/** @} */
-
-/** @defgroup aliases Common material aliases
- *  @{
- */
-#define DEDX_WATER DEDX_WATER_LIQUID
-#define DEDX_AIR DEDX_AIR_DRY_NEAR_SEA_LEVEL
-#define DEDX_PMMA DEDX_LUCITE_PERSPEX_PMMA
-#define DEDX_PERSPEX DEDX_LUCITE_PERSPEX_PMMA
-#define DEDX_LUCITE DEDX_LUCITE_PERSPEX_PMMA
-#define DEDX_TEFLON DEDX_POLYTETRAFLUOROETHYLENE
-#define DEDX_CONCRETE DEDX_CONCRETE_PORTLAND
-#define DEDX_CAESIUM DEDX_CESIUM
-/** @} */
-
 /** @brief Translate a numeric error code to a human-readable string.
  *  @param[out] err_str  Buffer to receive the error description (caller-allocated).
  *  @param[in]  err      Error code returned by a libdedx function.
@@ -483,32 +160,7 @@ float dedx_get_min_energy(int program, int ion);
  */
 float dedx_get_max_energy(int program, int ion);
 
-/** @cond INTERNAL */
-typedef struct {
-    int cache;
-    int hits;
-    int miss;
-} _dedx_lookup_accelerator;
-
-typedef struct {
-    float a;
-    float b;
-    float c;
-    float d;
-    float x;
-} _dedx_spline_base;
-
-typedef struct {
-    _dedx_spline_base base[_DEDX_MAXELEMENTS];
-    int n;
-    int prog;
-    int target;
-    int ion;
-    int datapoints;
-    _dedx_lookup_accelerator acc;
-} _dedx_lookup_data;
-
-/** @endcond */
+typedef struct _dedx_lookup_data _dedx_lookup_data;
 
 /**
  * @brief Workspace holding preloaded stopping power datasets.
@@ -585,8 +237,9 @@ typedef struct {
  *  @param[in]     ws      Workspace to load into.
  *  @param[in,out] config  Configuration to load; updated with resolved values.
  *  @param[out]    err     Error code; 0 on success.
+ *  @return 0 on success, -1 on failure.
  */
-void dedx_load_config(dedx_workspace *ws, dedx_config *config, int *err);
+int dedx_load_config(dedx_workspace *ws, dedx_config *config, int *err);
 
 /** @brief Evaluate the stopping power at a given energy.
  *
@@ -621,4 +274,4 @@ float dedx_get_simple_stp(int ion, int target, float energy, int *err);
  */
 void dedx_free_config(dedx_config *config, int *err);
 
-#endif // DEDX_H_INCLUDED
+#endif // DEDX_H
diff --git a/libdedx/dedx_atima.c b/libdedx/dedx_atima.c
deleted file mode 100644
index 47693c5..0000000
--- a/libdedx/dedx_atima.c
+++ /dev/null
@@ -1,192 +0,0 @@
-#include "dedx_atima.h"
-float _dedx_bethek(float energy, int PZ, int PA, int TZ, float TA, float pot, float rho);
-float _dedx_dedxela(float energy, int PZ, int PA, int TZ, float TA);
-float _dedx_sezi(float energy, int PZ, int PA, int TZ, float TA);
-float _dedx_rpstop(int TZ, float energy);
-float _dedx_bf4(int TZ, float eta, float zpeff);
-
-float _dedx_calculate_atima_energy(float energy, dedx_config *config) {
-    int *err = 0;
-    float PZ, PA, TZ, TA, rho, pot;
-    float dedx = 0;
-    float factor;
-
-    PZ = _dedx_get_atom_charge(config->ion, err);
-    PA = _dedx_get_atom_mass(config->ion, err);
-    TZ = _dedx_get_atom_charge(config->target, err);
-    TA = _dedx_get_atom_mass(config->target, err);
-    float e_tot = PA * energy * 1000;
-
-    if (energy <= 10) {
-        dedx = _dedx_sezi(e_tot, PZ, PA, TZ, TA) + _dedx_dedxela(e_tot, PZ, PA, TZ, TA);
-    } else if (energy <= 30) {
-        factor = 0.05 * (energy - 10);
-        dedx = (1 - factor) * _dedx_sezi(e_tot, PZ, PA, TZ, TA) + factor * _dedx_bethek(e_tot, PZ, PA, TZ, TA, pot, rho)
-               + _dedx_dedxela(e_tot, PZ, PA, TZ, TA);
-    } else {
-        dedx = _dedx_bethek(e_tot, PZ, PA, TZ, TA, pot, rho) + _dedx_dedxela(e_tot, PZ, PA, TZ, TA);
-    }
-    return dedx;
-}
-
-float _dedx_bethek(float energy, int PZ, int PA, int TZ, float TA, float pot, float rho) {
-    float e = energy / PA;
-    float e1 = energy / (1000 * PA);
-    float gamma, bsq, zeta, eta, zpeff, f1, f2, f4, f6;
-    float c, beta;
-
-    gamma = 1 + e1 / DEDX_AMU;
-    bsq = 1 - 1 / pow(gamma, 2);
-    beta = sqrt(bsq);
-    zeta = 1 - exp(-130 * beta / (pow(PZ, 2 / 3)));
-    zpeff = zeta * PZ;
-    eta = beta * gamma;
-
-    f1 = 4e-5 * DEDX_PI * DEDX_ECHARGE2 * DEDX_ECHARGE2 * DEDX_N_AVO / DEDX_EMASS * pow(zpeff, 2) * TZ / (bsq * TA);
-    f2 = log(2 * DEDX_EMASS * 1e6 * bsq / pot);
-    if (eta >= 0.13) {
-        c = (.422377 * pow(eta, -2) + .0304043 * pow(eta, -4) - .00038106 * pow(eta, -6)) * 1e-6 * pow(pot, 2)
-            + (3.858019 * pow(eta, -2) - 0.1667989 * (pow(eta, -4)) + 0.00157955 * (pow(eta, -6))) * 1.0e-9
-                  * pow(pot, 3);
-        f2 = f2 - c / TZ;
-    }
-    f4 = _dedx_bf4(TZ, eta, zpeff);
-    f6 = 2 * log(gamma) - bsq;
-    return 0;
-}
-
-float _dedx_bf4(int TZ, float eta, float zpeff) {
-    return 0;
-}
-
-float _dedx_dedxela(float energy, int PZ, int PA, int TZ, float TA) {
-    float sn, a;
-    float epsil = 32.53 * TA * energy / (PZ * TZ * (PA + TA) * (pow(PZ, 0.23) + pow(TZ, .23)));
-    if (epsil >= 30) {
-        sn = log(epsil) / (2 * epsil);
-    } else {
-        a = .01321 * pow(epsil, 0.21226) + .19593 * pow(epsil, 0.5);
-        sn = 0.5 * log(1 + 1.1383 * epsil) / (epsil + a);
-    }
-    sn = sn * PZ * TZ * PA * 8.462 / ((PA + TA) * pow(PZ, 0.23) + pow(TZ, 0.23));
-    sn = sn * .1 * DEDX_N_AVO / TA;
-    return sn;
-}
-
-float _dedx_sezi(float energy, int PZ, int PA, int TZ, float TA) {
-    float dedx = 0;
-    ;
-    int err = 0;
-    float *ion_info;
-    ion_info = _dedx_get_atima_data(PZ, &err);
-    float *target_info;
-    target_info = _dedx_get_atima_data(TZ, &err);
-    float v_fermi = target_info[1];
-    float lfctr = ion_info[2];
-    float e = energy / PA;
-    if (PZ == 1) {
-        dedx = _dedx_rpstop(TZ, e) * 0.60222 / TA;
-    } else if (PZ == 2) {
-        float E = e;
-        float he0 = 1;
-        float B, A, heh;
-        if (e < he0)
-            E = he0;
-        B = log(E);
-        A = 0.2865 + 0.1266 * B - .001429 * B * B + 0.02402 * B * B * B - 0.01135 * B * B * B * B
-            + 0.001475 * B * B * B * B * B;
-        if (A < 30)
-            heh = 1 - exp(-1 * A);
-        else
-            heh = 1 - exp(-30);
-        if (log(E) < 0)
-            A = (1.0 + (0.007 + 0.00005 * TZ) * exp(-pow(7.6, 2)));
-        else
-            A = (1.0 + (0.007 + 0.00005 * TZ) * exp(-pow(7.6 - log(E), 2)));
-        heh = heh * A * A;
-        dedx = _dedx_rpstop(TZ, e) * heh * 4;
-        if (e <= he0)
-            dedx = dedx * sqrt(e / he0);
-        dedx = dedx * 0.602222 / TA;
-    } else if (PZ > 2) {
-        float vr, yr, help, power, sp;
-        float a, q, b, l, l0, l1, q1, eee;
-        float zeta;
-        float yrmin = 0.130;
-        float vrmin = 1.0;
-        float vmin, help1, help2, help3, help4, help5;
-        float v = sqrt(e / 25.0) / v_fermi;
-        float v_2 = v * v;
-        if (v > 1.0) {
-            vr = v * v_fermi * (1 + 1 / (5 * v_2));
-        } else {
-            vr = 3 * v_fermi / 4 * (1 + (2 / 3 * v_2) - (v_2 * v_2 / 15));
-        }
-        yr = _dedx_max(vr / pow(PZ, 0.67), _dedx_max(vrmin / pow(PZ, 0.67), yrmin));
-        help = vr / pow(PZ, 0.6667);
-        yr = _dedx_max(yrmin, help);
-        help = vrmin / (pow(PZ, 0.6667));
-        yr = _dedx_max(yr, help);
-        a = -0.803 * pow(yr, 0.3) + 1.3167 * pow(yr, 0.6) + 00.38157 * yr + 0.008983 * yr * yr;
-        q = _dedx_min(1, _dedx_max(0, (1 - exp(-1 * _dedx_min(a, 50)))));
-        b = (_dedx_min(0.43, _dedx_max(0.32, 0.12 + 0.25 * PZ))) / pow(PZ, 0.3333);
-        l0 = (0.8 - q * _dedx_min(1.2, 0.6 + PZ / 30)) / pow(PZ, 0.3333);
-        if (q < 0.2) {
-            l1 = 0;
-        } else if (q < _dedx_max(0, 0.9 - 0.025 * PZ)) {
-            q1 = 0.2;
-            l1 = b * (q - 0.2) / fabs(_dedx_max(0, 0.9 - 0.025 * PZ) - 0.2000001);
-        } else if (q < _dedx_max(0, 1 - 0.025 * _dedx_min(16, PZ))) {
-            l1 = b;
-        } else {
-            l1 = b * (1 - q) / (0.025 * _dedx_min(16, PZ));
-        }
-        l = _dedx_max(l1, l0 * lfctr);
-        zeta = q + (1 / (2 * v_fermi * v_fermi)) * (1 - q) * log(1 + pow((4 * l * v_fermi / 1.919), 2));
-        a = -pow(7.6 - _dedx_max(0, log(e)), 2);
-        zeta = zeta * (1 + 1 / pow(PZ, 2)) * (0.18 + 0.0015 * TZ) * exp(a);
-        if (yr <= _dedx_max(yrmin, vrmin / pow(PZ, 0.6667))) {
-            vmin = 0.5 * (vrmin + sqrt(_dedx_max(0, vrmin * vrmin - 0.8 * v_fermi * v_fermi)));
-            eee = 25 * vmin * vmin;
-            power = 0.5;
-            if (TZ == 6 || ((TZ == 14 || TZ == 32) && PZ <= 19)) {
-                power = 0.35;
-            }
-            sp = _dedx_rpstop(TZ, e);
-            help1 = zeta * PZ;
-            help2 = help1 * help1;
-            help3 = e / eee;
-            help4 = pow(help3, power);
-            help5 = help2 * help4;
-            dedx = sp * help5;
-            dedx = dedx * 0.60222 / TA;
-        } else {
-            dedx = _dedx_rpstop(TZ, e) * (pow(zeta * PZ, 2) * 0.60222) / TA;
-        }
-    }
-
-    return dedx;
-}
-
-float _dedx_rpstop(int TZ, float energy) {
-    float sl, sh, sp;
-    float pe0 = 25;
-    float velpwr;
-    float E = energy;
-    const float *coef;
-    coef = dedx_pcoef[TZ];
-    if (E < pe0)
-        E = pe0;
-    sl = coef[0] * pow(E, coef[1]) + coef[2] * pow(E, coef[3]);
-    sh = coef[4] / (pow(E, coef[5])) * (log(coef[6] / E + coef[7] * E));
-    sp = sl * sh / (sl + sh);
-    if (E <= pe0) {
-        if (TZ > 6) {
-            velpwr = 0.45;
-        } else {
-            velpwr = 0.25;
-        }
-        sp = sp * pow((E / pe0), velpwr);
-    }
-    return sp;
-}
diff --git a/libdedx/dedx_atima.h b/libdedx/dedx_atima.h
deleted file mode 100644
index a847656..0000000
--- a/libdedx/dedx_atima.h
+++ /dev/null
@@ -1,112 +0,0 @@
-#ifndef DEDX_ATIMA_H_INCLUDED
-#define DEDX_ATIMA_H_INCLUDED
-#include <math.h>
-#include <stdio.h>
-#include <stdlib.h>
-#include <string.h>
-
-#include "dedx.h"
-#include "dedx_config.h"
-#include "dedx_const.h"
-#include "dedx_file.h"
-#include "dedx_file_access.h"
-#include "dedx_periodic_table.h"
-#include "dedx_split.h"
-#include "dedx_stopping_data.h"
-#include "tools/dedx_math.h"
-float _dedx_calculate_atima_energy(float energy, dedx_config *config);
-
-const float dedx_pcoef[100][8] = {{-1},
-                                  {.0091827, .0053496, .69741, .48493, 316.07, 1.0143, 9329.3, .0539890},
-                                  {.11393, .0051984, 1.0822, .39252, 1081.0, 1.0645, 4068.5, .0176990},
-                                  {.85837, .0050147, 1.6044, .38844, 1337.3, 1.047, 2659.2, .01898},
-                                  {.8781, .0051049, 5.4232, .2032, 1200.6, 1.0211, 1401.8, .0385290},
-                                  {1.4608, .0048836, 2.338, .44249, 1801.3, 1.0352, 1784.1, .02024},
-                                  {3.2579, .0049148, 2.7156, .36473, 2092.2, 1.0291, 2643.6, .0182370},
-                                  {.59674, .0050837, 4.2073, .30612, 2394.2, 1.0255, 4892.1, .0160060},
-                                  {.75253, .0050314, 4.0824, .30067, 2455.8, 1.0181, 5069.7, .0174260},
-                                  {1.226, .0051385, 3.2246, .32703, 2525.0, 1.0142, 7563.6, .0194690},
-                                  {1.0332, .0051645, 3.004, .33889, 2338.6, .99997, 6991.2, .0217990},
-                                  {6.0972, .0044292, 3.1929, .45763, 1363.3, .95182, 2380.6, .0818350},
-                                  {14.013, .0043646, 2.2641, .36326, 2187.4, .99098, 6264.8, .0462},
-                                  {.039001, .0045415, 5.5463, .39562, 1589.2, .95316, 816.16, .0474840},
-                                  {2.072, .0044516, 3.5585, .53933, 1515.2, .93161, 1790.3, .0351980},
-                                  {17.575, .0038346, .078694, 1.2388, 2806.0, .97284, 1037.6, .0128790},
-                                  {16.126, .0038315, .054164, 1.3104, 2813.3, .96587, 1251.4, .0118470},
-                                  {3.217, .0044579, 3.6696, .5091, 2734.6, .96253, 2187.5, .0169070},
-                                  {2.0379, .0044775, 3.0743, .54773, 3505.0, .97575, 1714.00, .0117010},
-                                  {.74171, .0043051, 1.1515, .95083, 917.21, .8782, 389.93, .18926},
-                                  {9.1316, .0043809, 5.4611, .31327, 3891.8, .97933, 6267.9, .0151960},
-                                  {7.2247, .0043718, 6.1017, .37511, 2829.2, .95218, 6376.1, .0203980},
-                                  {.147, .0048456, 6.3485, .41057, 2164.1, .94028, 5292.6, .0502630},
-                                  {5.0611, .0039867, 2.6174, .57957, 2218.9, .92361, 6323.00, .0256690},
-                                  {.53267, .0042968, .39005, 1.2725, 1872.7, .90776, 64.166, .0301070},
-                                  {.47697, .0043038, .31452, 1.3289, 1920.5, .90649, 45.576, .0274690},
-                                  {.027426, .0035443, .031563, 2.1755, 1919.5, .90099, 23.902, .0253630},
-                                  {.16383, .0043042, .073454, 1.8592, 1918.4, .89678, 27.61, .0231840},
-                                  {4.2562, .0043737, 1.5606, .72067, 1546.8, .87958, 302.02, .0409440},
-                                  {2.3508, .0043237, 2.882, .50113, 1837.7, .89992, 2377.00, .04965},
-                                  {3.1095, .0038455, .11477, 1.5037, 2184.7, .89309, 67.306, .0165880},
-                                  {15.322, .0040306, .65391, .67668, 3001.7, .92484, 3344.2, .0163660},
-                                  {3.6932, .0044813, 8.608, .27638, 2982.7, .9276, 3166.6, .0308740},
-                                  {7.1373, .0043134, 9.4247, .27937, 2725.8, .91597, 3166.1, .0250080},
-                                  {4.8979, .0042937, 3.7793, .50004, 2824.5, .91028, 1282.4, .0170610},
-                                  {1.3683, .0043024, 2.5679, .60822, 6907.8, .9817, 628.01, .0068055},
-                                  {1.8301, .0042983, 2.9057, .6038, 4744.6, .94722, 936.64, .0092242},
-                                  {.42056, .0041169, .01695, 2.3616, 2252.7, .89192, 39.752, .0277570},
-                                  {30.78, .0037736, .55813, .76816, 7113.2, .97697, 1604.4, .0065268},
-                                  {11.576, .0042119, 7.0244, .37764, 4713.5, .94264, 2493.2, .01127},
-                                  {6.2406, .0041916, 5.2701, .49453, 4234.6, .93232, 2063.9, .0118440},
-                                  {.33073, .0041243, 1.7246, 1.1062, 1930.2, .86907, 27.416, .0382080},
-                                  {.017747, .0041715, .14586, 1.7305, 1803.6, .86315, 29.669, .0321230},
-                                  {3.7229, .0041768, 4.6286, .56769, 1678.0, .86202, 3094.00, .06244},
-                                  {.13998, .0041329, .25573, 1.4241, 1919.3, .86326, 72.797, .0322350},
-                                  {.2859, .0041386, .31301, 1.3424, 1954.8, .86175, 115.18, .0293420},
-                                  {.76, .0042179, 3.386, .76285, 1867.4, .85805, 69.994, .0364480},
-                                  {6.3957, .0041935, 5.4689, .41378, 1712.6, .85397, 18493.00, .0564710},
-                                  {3.4717, .0041344, 3.2337, .63788, 1116.4, .81959, 4766.0, .1179},
-                                  {2.5265, .0042282, 4.532, .53562, 1030.8, .81652, 16252.0, .19722},
-                                  {7.3683, .0041007, 4.6791, .51428, 1160.0, .82454, 17965.0, .13316},
-                                  {7.7197, .004388, 3.242, .68434, 1428.1, .83398, 1786.7, .0665120},
-                                  {16.78, .0041918, 9.3198, .29568, 3370.9, .90289, 7431.7, .02616},
-                                  {4.2132, .0042098, 4.6753, .57945, 3503.9, .89261, 1468.9, .0143590},
-                                  {4.0818, .004214, 4.4425, .58393, 3945.3, .90281, 1340.5, .0134140},
-                                  {.18517, .0036215, .00058788, 3.5315, 2931.3, .88936, 26.18, .0263930},
-                                  {4.8248, .0041458, 6.0934, .57026, 2300.1, .86359, 2980.7, .0386790},
-                                  {.49857, .0041054, 1.9775, .95877, 786.55, .78509, 806.6, .40882},
-                                  {3.2754, .0042177, 5.768, .54054, 6631.3, .94282, 744.07, .0083026},
-                                  {2.9978, .0040901, 4.5299, .62025, 2161.2, .85669, 1268.6, .0430310},
-                                  {2.8701, .004096, 4.2568, .6138, 2130.4, .85235, 1704.1, .0393850},
-                                  {10.853, .0041149, 5.8907, .46834, 2857.2, .8755, 3654.2, .0299550},
-                                  {3.6407, .0041782, 4.8742, .57861, 1267.7, .82211, 3508.2, .24174},
-                                  {17.645, .0040992, 6.5855, .32734, 3931.3, .90754, 5156.7, .0362780},
-                                  {7.5309, .0040814, 4.9389, .50679, 2519.7, .85819, 3314.6, .0305140},
-                                  {5.4742, .0040829, 4.897, .51113, 2340.1, .85296, 2342.7, .0356620},
-                                  {4.2661, .0040667, 4.5032, .55257, 2076.4, .84151, 1666.6, .0408010},
-                                  {6.8313, .0040486, 4.3987, .51675, 2003.0, .83437, 1410.4, .03478},
-                                  {1.2707, .0040553, 4.6295, .57428, 1626.3, .81858, 995.68, .0553190},
-                                  {5.7561, .0040491, 4.357, .52496, 2207.3, .83796, 1579.5, .0271650},
-                                  {14.127, .0040596, 5.8304, .37755, 3645.9, .87823, 3411.8, .0163920},
-                                  {6.6948, .0040603, 4.9361, .47961, 2719.0, .85249, 1885.8, .0197130},
-                                  {3.0619, .0040511, 3.5803, .59082, 2346.1, .83713, 1222.0, .0200720},
-                                  {10.811, .0033008, 1.3767, .76512, 2003.7, .82269, 1110.6, .0249580},
-                                  {2.7101, .0040961, 1.2289, .98598, 1232.4, .79066, 155.42, .0472940},
-                                  {.52345, .0040244, 1.4038, .8551, 1461.4, .79677, 503.34, .0367890},
-                                  {.4616, .0040203, 1.3014, .87043, 1473.5, .79687, 443.09, .0363010},
-                                  {.97814, .0040374, 2.0127, .7225, 1890.8, .81747, 930.7, .02769},
-                                  {3.2086, .004051, 3.6658, .53618, 3091.2, .85602, 1508.1, .0154010},
-                                  {2.0035, .0040431, 7.4882, .3561, 4464.3, .88836, 3966.5, .0128390},
-                                  {15.43, .0039432, 1.1237, .70703, 4595.7, .88437, 1576.5, .0088534},
-                                  {3.1512, .0040524, 4.0996, .5425, 3246.3, .85772, 1691.8, .0150580},
-                                  {7.1896, .0040588, 8.6927, .35842, 4760.6, .88833, 2888.3, .0110290},
-                                  {9.3209, .004054, 11.543, .32027, 4866.2, .89124, 3213.4, .0119350},
-                                  {29.242, .0036195, .16864, 1.1226, 5688.0, .89812, 1033.3, .0071303},
-                                  {1.8522, .0039973, 3.1556, .65096, 3755.0, .86383, 1602.0, .0120420},
-                                  {3.222, .0040041, 5.9024, .52678, 4040.2, .86804, 1658.4, .0117470},
-                                  {9.3412, .0039661, 7.921, .42977, 5180.9, .88773, 2173.2, .0092007},
-                                  {36.183, .0036003, .58341, .86747, 6990.2, .91082, 1417.1, .0062187},
-                                  {5.9284, .0039695, 6.4082, .52122, 4619.5, .88083, 2323.5, .0116270},
-                                  {5.2454, .0039744, 6.7969, .48542, 4586.3, .87794, 2481.5, .0112820},
-                                  {33.702, .0036901, .47257, .89235, 5295.7, .8893, 2053.3, .0091908},
-                                  {2.7589, .0039806, 3.2092, .66122, 2505.4, .82863, 2065.1, .0228160}};
-#endif
diff --git a/libdedx/dedx_bethe.c b/libdedx/dedx_bethe.c
index 5d71573..91a2c09 100644
--- a/libdedx/dedx_bethe.c
+++ b/libdedx/dedx_bethe.c
@@ -16,12 +16,17 @@
 #include "dedx_bethe.h"
 
 #include <math.h>
+#include <stdlib.h>
 
-void _dedx_gold_section(_dedx_bethe_struct bet, _dedx_gold_struct *gold, int *err);
-float _dedx_m(float PT, _dedx_bethe_struct bet, int *err);
-float _dedx_mm(float PT, _dedx_bethe_struct bet, _dedx_gold_struct gold, int *err);
+static void gold_section(_dedx_bethe_struct bet, _dedx_gold_struct *gold, int *err);
 
-float _dedx_calculate_bethe_energy(
+/* Raw Bethe-model evaluator without the stitched low-energy extension. */
+static float evaluate_bethe_model(float PT, _dedx_bethe_struct bet, int *err);
+
+/* Full evaluator with Lindhard-Scharff behavior below the sewn transition. */
+static float evaluate_bethe_model_LEext(float PT, _dedx_bethe_struct bet, _dedx_gold_struct gold, int *err);
+
+float dedx_internal_calculate_bethe_energy(
     _dedx_bethe_coll_struct *ws, float energy, float PZ, float PA, float TZ, float TA, float rho, float Io_Pot) {
 
     float dedx;
@@ -31,7 +36,7 @@ float _dedx_calculate_bethe_energy(
     if (gold != NULL && bet != NULL
         && (bet->PZ0 == PZ && bet->PA0 == PA && bet->TZ0 == TZ && bet->TA0 == TA && bet->rho == rho
             && bet->potentiale == Io_Pot)) {
-        dedx = _dedx_mm(energy * PA, *bet, *gold, &err);
+        dedx = evaluate_bethe_model_LEext(energy * PA, *bet, *gold, &err);
         return dedx;
     } else {
         if (ws->gold != NULL)
@@ -59,14 +64,14 @@ float _dedx_calculate_bethe_energy(
     bet->TZ0 = TZ;
     bet->PA0 = PA;
 
-    _dedx_gold_section(*bet, gold, &err);
+    gold_section(*bet, gold, &err);
 
-    dedx = _dedx_mm(energy * PA, *bet, *gold, &err);
+    dedx = evaluate_bethe_model_LEext(energy * PA, *bet, *gold, &err);
 
     return dedx;
 }
 
-float _dedx_mm(float PT, _dedx_bethe_struct bet, _dedx_gold_struct gold, int *err) {
+static float evaluate_bethe_model_LEext(float PT, _dedx_bethe_struct bet, _dedx_gold_struct gold, int *err) {
     double T = PT;
     float dedx;
     double mass = 940 * bet.PA0;
@@ -135,19 +140,18 @@ float _dedx_mm(float PT, _dedx_bethe_struct bet, _dedx_gold_struct gold, int *er
     return dedx;
 }
 
-void _dedx_gold_section(_dedx_bethe_struct bet, _dedx_gold_struct *gold, int *err) {
+static void gold_section(_dedx_bethe_struct bet, _dedx_gold_struct *gold, int *err) {
 
     double e_min = gold->e_min;
     double e_max = gold->e_max;
-    double e_zero = gold->e_zero;
 
     double rla = e_min + 0.3819661 * (e_max - e_min);
     double rmu = e_min + 0.6180339 * (e_max - e_min);
-    double tla = _dedx_m(rla * bet.PA0, bet, err);
-    double tmu = _dedx_m(rmu * bet.PA0, bet, err);
+    double tla = evaluate_bethe_model(rla * bet.PA0, bet, err);
+    double tmu = evaluate_bethe_model(rmu * bet.PA0, bet, err);
     while (1) {
         if ((e_max - e_min) < gold->epsilon) {
-            e_zero = (e_max + e_min) * 0.5;
+            gold->e_zero = (e_max + e_min) * 0.5;
             break;
         } else {
             if (tla > 0.0) {
@@ -156,30 +160,30 @@ void _dedx_gold_section(_dedx_bethe_struct bet, _dedx_gold_struct *gold, int *er
                 tmu = tla;
                 rla = e_min + 0.3819661 * (e_max - e_min);
 
-                tla = _dedx_m(rla * bet.PA0, bet, err);
+                tla = evaluate_bethe_model(rla * bet.PA0, bet, err);
             } else {
 
                 e_min = rla;
                 rla = rmu;
                 tla = tmu;
                 rmu = e_min + 0.6180339 * (e_max - e_min);
-                tmu = _dedx_m(rmu * bet.PA0, bet, err);
+                tmu = evaluate_bethe_model(rmu * bet.PA0, bet, err);
             }
         }
     }
 
-    e_min = e_zero;
+    e_min = gold->e_zero;
     e_max = gold->e_max;
 
     rla = e_min + 0.3819661 * (e_max - e_min);
     rmu = e_min + 0.6180339 * (e_max - e_min);
-    tla = _dedx_m(rla * bet.PA0, bet, err);
-    tmu = _dedx_m(rmu * bet.PA0, bet, err);
+    tla = evaluate_bethe_model(rla * bet.PA0, bet, err);
+    tmu = evaluate_bethe_model(rmu * bet.PA0, bet, err);
 
     while (1) {
         if ((e_max - e_min) < gold->epsilon) {
             gold->e_extr = (e_min + e_max) * 0.5;
-            gold->f_extr = _dedx_m(gold->e_extr * bet.PA0, bet, err);
+            gold->f_extr = evaluate_bethe_model(gold->e_extr * bet.PA0, bet, err);
             break;
         } else {
             if (tla > tmu) {
@@ -187,13 +191,13 @@ void _dedx_gold_section(_dedx_bethe_struct bet, _dedx_gold_struct *gold, int *er
                 rmu = rla;
                 tmu = tla;
                 rla = e_min + 0.3819661 * (e_max - e_min);
-                tla = _dedx_m(rla * bet.PA0, bet, err);
+                tla = evaluate_bethe_model(rla * bet.PA0, bet, err);
             } else {
                 e_min = rla;
                 rla = rmu;
                 tla = tmu;
                 rmu = e_min + 0.6180339 * (e_max - e_min);
-                tmu = _dedx_m(rmu * bet.PA0, bet, err);
+                tmu = evaluate_bethe_model(rmu * bet.PA0, bet, err);
             }
         }
     }
@@ -202,14 +206,18 @@ void _dedx_gold_section(_dedx_bethe_struct bet, _dedx_gold_struct *gold, int *er
 
     rla = e_min + 0.3819661 * (e_max - e_min);
     rmu = e_min + 0.6180339 * (e_max - e_min);
-    tla = (_dedx_m((rla + gold->h) * bet.PA0, bet, err) - _dedx_m((rla - gold->h) * bet.PA0, bet, err)) / gold->h
-          - _dedx_m(rla * bet.PA0, bet, err) / rla;
-    tmu = (_dedx_m((rmu + gold->h) * bet.PA0, bet, err) - _dedx_m((rmu - gold->h) * bet.PA0, bet, err)) / gold->h
-          - _dedx_m(rmu * bet.PA0, bet, err) / rmu;
+    tla = (evaluate_bethe_model((rla + gold->h) * bet.PA0, bet, err)
+           - evaluate_bethe_model((rla - gold->h) * bet.PA0, bet, err))
+              / gold->h
+          - evaluate_bethe_model(rla * bet.PA0, bet, err) / rla;
+    tmu = (evaluate_bethe_model((rmu + gold->h) * bet.PA0, bet, err)
+           - evaluate_bethe_model((rmu - gold->h) * bet.PA0, bet, err))
+              / gold->h
+          - evaluate_bethe_model(rmu * bet.PA0, bet, err) / rmu;
     while (1) {
         if ((e_max - e_min) < gold->epsilon) {
             gold->e_sewn = (e_min + e_max) * 0.5;
-            gold->f_sewn = _dedx_m(gold->e_sewn * bet.PA0, bet, err);
+            gold->f_sewn = evaluate_bethe_model(gold->e_sewn * bet.PA0, bet, err);
             break;
         } else {
             if (tla <= 0) {
@@ -217,24 +225,26 @@ void _dedx_gold_section(_dedx_bethe_struct bet, _dedx_gold_struct *gold, int *er
                 rmu = rla;
                 tmu = tla;
                 rla = e_min + 0.3819661 * (e_max - e_min);
-                tla = (_dedx_m((rla + gold->h) * bet.PA0, bet, err) - _dedx_m((rla - gold->h) * bet.PA0, bet, err))
+                tla = (evaluate_bethe_model((rla + gold->h) * bet.PA0, bet, err)
+                       - evaluate_bethe_model((rla - gold->h) * bet.PA0, bet, err))
                           / gold->h
-                      - _dedx_m(rla * bet.PA0, bet, err) / rla;
+                      - evaluate_bethe_model(rla * bet.PA0, bet, err) / rla;
 
             } else {
                 e_min = rla;
                 rla = rmu;
                 tla = tmu;
                 rmu = e_min + 0.6180339 * (e_max - e_min);
-                tmu = (_dedx_m((rmu + gold->h) * bet.PA0, bet, err) - _dedx_m((rmu - gold->h) * bet.PA0, bet, err))
+                tmu = (evaluate_bethe_model((rmu + gold->h) * bet.PA0, bet, err)
+                       - evaluate_bethe_model((rmu - gold->h) * bet.PA0, bet, err))
                           / gold->h
-                      - _dedx_m(rmu * bet.PA0, bet, err) / rmu;
+                      - evaluate_bethe_model(rmu * bet.PA0, bet, err) / rmu;
             }
         }
     }
 }
 
-float _dedx_m(float PT, _dedx_bethe_struct bet, int *err) {
+static float evaluate_bethe_model(float PT, _dedx_bethe_struct bet, int *err) {
     double T = PT;
     double mass = 940 * bet.PA0;
 
diff --git a/libdedx/dedx_bethe.h b/libdedx/dedx_bethe.h
index bf8b03c..594a8b7 100644
--- a/libdedx/dedx_bethe.h
+++ b/libdedx/dedx_bethe.h
@@ -1,39 +1,20 @@
-#ifndef DEDX_BETHE_H_INCLUDED
-#define DEDX_BETHE_H_INCLUDED
+#ifndef DEDX_BETHE_H
+#define DEDX_BETHE_H
 
-#include "dedx.h"
-#include "dedx_file_access.h"
+#include "dedx_bethe_struct.h"
 
-// #include "dedx_gold_struct.h"
-// #include "dedx_bethe_struct.h"
-
-typedef struct {
-    double TZ0;
-    double TA0;
-    double potentiale;
-    double rho;
-    double PZ0;
-    double PA0;
-} _dedx_bethe_struct;
-
-typedef struct {
-    double e_min;
-    double e_max;
-    double epsilon;
-    double h;
-    double e_zero;
-    double e_extr;
-    double f_extr;
-    double e_sewn;
-    double f_sewn;
-} _dedx_gold_struct;
-
-typedef struct {
-    _dedx_bethe_struct *bet;
-    _dedx_gold_struct *gold;
-} _dedx_bethe_coll_struct;
-
-float _dedx_calculate_bethe_energy(
+/** @brief Evaluate Bethe stopping power for one energy using cached work buffers.
+ *  @param[in,out] ws      Cached Bethe/Gauss-search workspace.
+ *  @param[in]     energy  Projectile energy in MeV/u.
+ *  @param[in]     PZ      Projectile charge.
+ *  @param[in]     PA      Projectile mass number.
+ *  @param[in]     TZ      Target atomic number.
+ *  @param[in]     TA      Target mass number.
+ *  @param[in]     rho     Target density in g/cm^3.
+ *  @param[in]     Io_Pot  Mean excitation potential in eV.
+ *  @return Stopping power in MeV cm^2/g.
+ */
+float dedx_internal_calculate_bethe_energy(
     _dedx_bethe_coll_struct *ws, float energy, float PZ, float PA, float TZ, float TA, float rho, float Io_Pot);
 
-#endif // DEDX_BETHE_H_INCLUDED
+#endif // DEDX_BETHE_H
diff --git a/libdedx/dedx_bethe_struct.h b/libdedx/dedx_bethe_struct.h
index 309f7c1..e34cb5f 100644
--- a/libdedx/dedx_bethe_struct.h
+++ b/libdedx/dedx_bethe_struct.h
@@ -1,13 +1,33 @@
-#ifndef DEDX_BETHE_STRUCT_H_INCLUDED
-#define DEDX_BETHE_STRUCT_H_INCLUDED
-/*typedef struct
-{
+#ifndef DEDX_BETHE_STRUCT_H
+#define DEDX_BETHE_STRUCT_H
+
+/** @brief Cached Bethe input parameters for one target/projectile combination. */
+typedef struct {
     double TZ0;
     double TA0;
     double potentiale;
     double rho;
     double PZ0;
     double PA0;
-} bethe_struct;*/
+} _dedx_bethe_struct;
+
+/** @brief Scratch values used while stitching the low-energy extension. */
+typedef struct {
+    double e_min;
+    double e_max;
+    double epsilon;
+    double h;
+    double e_zero;
+    double e_extr;
+    double f_extr;
+    double e_sewn;
+    double f_sewn;
+} _dedx_gold_struct;
+
+/** @brief Heap-owned Bethe work buffers cached across repeated evaluations. */
+typedef struct {
+    _dedx_bethe_struct *bet;
+    _dedx_gold_struct *gold;
+} _dedx_bethe_coll_struct;
 
-#endif // DEDX_BETHE_STRUCT_H_INCLUDED
+#endif // DEDX_BETHE_STRUCT_H
diff --git a/libdedx/dedx_calculate_charge.c b/libdedx/dedx_calculate_charge.c
deleted file mode 100644
index 95a658f..0000000
--- a/libdedx/dedx_calculate_charge.c
+++ /dev/null
@@ -1,28 +0,0 @@
-#include "dedx_calculate_charge.h"
-
-#include "dedx_periodic_table.h"
-
-float _dedx_get_avg_charge(float composition[][2], int length) {
-    int i = 0;
-    int err = 0;
-    float charge = 0.0;
-    float nucleon = 0.0;
-    for (i = 0; i < length; i++) {
-        nucleon += composition[i][1] / _dedx_get_atom_mass((int) composition[i][0], &err);
-    }
-    for (i = 0; i < length; i++) {
-        charge += composition[i][1] / _dedx_get_atom_mass((int) composition[i][0], &err) / nucleon
-                  * _dedx_get_atom_charge((int) composition[i][0], &err);
-    }
-    return charge;
-}
-
-float _dedx_get_weight_charge(float composition[][2], int length) {
-    int err = 0;
-    float charge = 0.0;
-    int i = 0;
-    for (i = 0; i < length; i++) {
-        charge += composition[i][1] * _dedx_get_atom_charge((int) composition[i][0], &err);
-    }
-    return charge;
-}
diff --git a/libdedx/dedx_calculate_charge.h b/libdedx/dedx_calculate_charge.h
deleted file mode 100644
index 4cbd71e..0000000
--- a/libdedx/dedx_calculate_charge.h
+++ /dev/null
@@ -1,7 +0,0 @@
-#ifndef DEDX_CALCULATE_CHARGE_H_INCLUDED
-#define DEDX_CALCULATE_CHARGE_H_INCLUDED
-
-float _dedx_get_avg_charge(float composition[][2], int length);
-float _dedx_get_weight_charge(float composition[][2], int length);
-
-#endif // DEDX_CALCULATE_CHARGE_H_INCLUDED
diff --git a/libdedx/dedx_const.h b/libdedx/dedx_const.h
index 11734ea..8673a86 100644
--- a/libdedx/dedx_const.h
+++ b/libdedx/dedx_const.h
@@ -1,10 +1,30 @@
-#ifndef DEDX_CONST_H_INCLUDED
-#define DEDX_CONST_H_INCLUDED
+#ifndef DEDX_CONST_H
+#define DEDX_CONST_H
+
+/**
+ * @file dedx_const.h
+ * @brief Shared physical constants used by libdedx internals.
+ *
+ * These are compile-time constants in the legacy units expected by the
+ * original stopping-power implementations bundled in libdedx.
+ */
+
+/** Avogadro constant in mol^-1. */
 #define DEDX_N_AVO 6.0221367
+
+/** Atomic mass unit in MeV/c^2. */
 #define DEDX_AMU 931.49432
+
+/** Pi. */
 #define DEDX_PI 3.14159265358979323846
+
+/** Electron rest mass in MeV/c^2. */
 #define DEDX_EMASS 0.51099906
+
+/** Elementary charge scaling used by the legacy formulas. */
 #define DEDX_ECHARGE 1.60217733
+
+/** e^2 / (4 pi epsilon_0) in eV*Angstrom units used by the Bethe code. */
 #define DEDX_ECHARGE2 14.399652
 
-#endif
+#endif // DEDX_CONST_H
diff --git a/libdedx/dedx_elements.h b/libdedx/dedx_elements.h
new file mode 100644
index 0000000..6f55de6
--- /dev/null
+++ b/libdedx/dedx_elements.h
@@ -0,0 +1,329 @@
+#ifndef DEDX_ELEMENTS_H
+#define DEDX_ELEMENTS_H
+
+/** Maximum number of tabulated energy points stored per dataset. */
+#define DEDX_MAX_ELEMENTS 150
+
+/**
+ * @defgroup ions_and_materials Ion and material identifiers
+ * @brief Identifiers for projectile ions and target materials.
+ *
+ * Elemental ions (Z=1–98) correspond to their atomic number.
+ * Compound materials follow sequentially after the elements.
+ * Both share a single enum to preserve the numeric ID mapping.
+ * @{
+ */
+enum {
+    /* --- Elemental ions (Z = atomic number) --- */
+    DEDX_HYDROGEN = 1,
+    DEDX_HELIUM,
+    DEDX_LITHIUM,
+    DEDX_BERYLLIUM,
+    DEDX_BORON,
+    DEDX_CARBON,
+    DEDX_GRAPHITE = 906,
+    DEDX_NITROGEN = 7,
+    DEDX_OXYGEN,
+    DEDX_FLUORINE,
+    DEDX_NEON,
+    DEDX_SODIUM,
+    DEDX_MAGNESIUM,
+    DEDX_ALUMINUM,
+    DEDX_SILICON,
+    DEDX_PHOSPHORUS,
+    DEDX_SULFUR,
+    DEDX_CHLORINE,
+    DEDX_ARGON,
+    DEDX_POTASSIUM,
+    DEDX_CALCIUM,
+    DEDX_SCANDIUM,
+    DEDX_TITANIUM,
+    DEDX_VANADIUM,
+    DEDX_CHROMIUM,
+    DEDX_MANGANESE,
+    DEDX_IRON,
+    DEDX_COBALT,
+    DEDX_NICKEL,
+    DEDX_COPPER,
+    DEDX_ZINC,
+    DEDX_GALLIUM,
+    DEDX_GERMANIUM,
+    DEDX_ARSENIC,
+    DEDX_SELENIUM,
+    DEDX_BROMINE,
+    DEDX_KRYPTON,
+    DEDX_RUBIDIUM,
+    DEDX_STRONTIUM,
+    DEDX_YTTRIUM,
+    DEDX_ZIRCONIUM,
+    DEDX_NIOBIUM,
+    DEDX_MOLYBDENUM,
+    DEDX_TECHNETIUM,
+    DEDX_RUTHENIUM,
+    DEDX_RHODIUM,
+    DEDX_PALLADIUM,
+    DEDX_SILVER,
+    DEDX_CADMIUM,
+    DEDX_INDIUM,
+    DEDX_TIN,
+    DEDX_ANTIMONY,
+    DEDX_TELLURIUM,
+    DEDX_IODINE,
+    DEDX_XENON,
+    DEDX_CESIUM,
+    DEDX_BARIUM,
+    DEDX_LANTHANUM,
+    DEDX_CERIUM,
+    DEDX_PRASEODYMIUM,
+    DEDX_NEODYMIUM,
+    DEDX_PROMETHIUM,
+    DEDX_SAMARIUM,
+    DEDX_EUROPIUM,
+    DEDX_GADOLINIUM,
+    DEDX_TERBIUM,
+    DEDX_DYSPROSIUM,
+    DEDX_HOLMIUM,
+    DEDX_ERBIUM,
+    DEDX_THULIUM,
+    DEDX_YTTERBIUM,
+    DEDX_LUTETIUM,
+    DEDX_HAFNIUM,
+    DEDX_TANTALUM,
+    DEDX_TUNGSTEN,
+    DEDX_RHENIUM,
+    DEDX_OSMIUM,
+    DEDX_IRIDIUM,
+    DEDX_PLATINUM,
+    DEDX_GOLD,
+    DEDX_MERCURY,
+    DEDX_THALLIUM,
+    DEDX_LEAD,
+    DEDX_BISMUTH,
+    DEDX_POLONIUM,
+    DEDX_ASTATINE,
+    DEDX_RADON,
+    DEDX_FRANCIUM,
+    DEDX_RADIUM,
+    DEDX_ACTINIUM,
+    DEDX_THORIUM,
+    DEDX_PROTACTINIUM,
+    DEDX_URANIUM,
+    DEDX_NEPTUNIUM,
+    DEDX_PLUTONIUM,
+    DEDX_AMERICIUM,
+    DEDX_CURIUM,
+    DEDX_BERKELIUM,
+    DEDX_CALIFORNIUM,
+
+    /* --- Compound and mixture targets (NIST material database) --- */
+    DEDX_A150_TISSUE_EQUIVALENT_PLASTIC,
+    DEDX_ACETONE,
+    DEDX_ACETYLENE,
+    DEDX_ADENINE,
+    DEDX_ADIPOSE_TISSUE_ICRP,
+    DEDX_AIR_DRY_NEAR_SEA_LEVEL,
+    DEDX_ALANINE,
+    DEDX_ALUMINUMOXIDE,
+    DEDX_AMBER,
+    DEDX_AMMONIA,
+    DEDX_ANILINE,
+    DEDX_ANTHRACENE,
+    DEDX_B100,
+    DEDX_BAKELITE,
+    DEDX_BARIUM_FLUORIDE,
+    DEDX_BARIUM_SULFATE,
+    DEDX_BENZENE,
+    DEDX_BERYLLIUM_OXIDE,
+    DEDX_BISMUTH_GERMANIUM_OXIDE,
+    DEDX_BLOOD_ICRP,
+    DEDX_BONE_COMPACT_ICRU,
+    DEDX_BONE_CORTICAL_ICRP,
+    DEDX_BORON_CARBIDE,
+    DEDX_BORON_OXIDE,
+    DEDX_BRAIN_ICRP,
+    DEDX_BUTANE,
+    DEDX_N_BUTYLALCOHOL,
+    DEDX_C552_AIR_EQUIVALENT_PLASTIC,
+    DEDX_CADMIUM_TELLURIDE,
+    DEDX_CADMIUM_TUNGSTATE,
+    DEDX_CALCIUM_CARBONATE,
+    DEDX_CALCIUM_FLUORIDE,
+    DEDX_CALCIUM_OXIDE,
+    DEDX_CALCIUM_SULFATE,
+    DEDX_CALCIUM_TUNGSTATE,
+    DEDX_CARBON_DIOXIDE,
+    DEDX_CARBON_TETRACHLORIDE,
+    DEDX_CELLULOSE_ACETATE_CELLOPHANE,
+    DEDX_CELLULOSE_ACETATE_BUTYRATE,
+    DEDX_CELLULOSE_NITRATE,
+    DEDX_CERIC_SULFATE_DOSIMETER_SOLUTION,
+    DEDX_CESIUM_FLUORIDE,
+    DEDX_CESIUM_IODIDE,
+    DEDX_CHLORO_BENZENE,
+    DEDX_CHLOROFORM,
+    DEDX_CONCRETE_PORTLAND,
+    DEDX_CYCLOHEXANE,
+    DEDX_DICHLOROBENZENE,
+    DEDX_DICHLORODIETHYL_ETHER,
+    DEDX_DICHLOROETHANE,
+    DEDX_DIETHYLETHER,
+    DEDX_N_N_DIMETHYL_FORMAMIDE,
+    DEDX_DIMETHYL_SULFOXIDE,
+    DEDX_ETHANE,
+    DEDX_ETHYL_ALCOHOL,
+    DEDX_ETHYL_CELLULOSE,
+    DEDX_ETHYLENE,
+    DEDX_EYE_LENS_ICRP,
+    DEDX_FERRIC_OXIDE,
+    DEDX_FERRO_BORIDE,
+    DEDX_FERROUS_OXIDE,
+    DEDX_FERROUS_SULFATE_DOSIMETER_SOLUTION,
+    DEDX_FREON_12,
+    DEDX_FREON_12B2,
+    DEDX_FREON_13,
+    DEDX_FREON_13B1,
+    DEDX_FREON_13I1,
+    DEDX_GADOLINIUM_OXYSULFIDE,
+    DEDX_GALLIUM_ARSENIDE,
+    DEDX_GEL_IN_PHOTOGRAPHIC_EMULSION,
+    DEDX_GLASS_PYREX,
+    DEDX_GLASS_LEAD,
+    DEDX_GLASS_PLATE, /* 169,170,171 */
+    DEDX_GLUCOSE,
+    DEDX_GLUTAMINE,
+    DEDX_GLYCEROL,
+    DEDX_GUANINE,
+    DEDX_GYPSUM_PLASTER_OF_PARIS,
+    DEDX_N_HEPTANE,
+    DEDX_N_HEXANE,
+    DEDX_KAPTON_POLYIMIDE_FILM,
+    DEDX_LANTHANUM_OXYBROMIDE,
+    DEDX_LANTHANUM_OXYSULFIDE,
+    DEDX_LEAD_OXIDE,
+    DEDX_LITHIUM_AMIDE,
+    DEDX_LITHIUM_CARBONATE,
+    DEDX_LITHIUM_FLUORIDE,
+    DEDX_LITHIUM_HYDRIDE,
+    DEDX_LITHIUM_IODIDE,
+    DEDX_LITHIUM_OXIDE,
+    DEDX_LITHIUM_TETRABORATE,
+    DEDX_LUNG_ICRP,
+    DEDX_M3_WAX,
+    DEDX_MAGNESIUM_CARBONATE,
+    DEDX_MAGNESIUM_FLUORIDE,
+    DEDX_MAGNESIUM_OXIDE,
+    DEDX_MAGNESIUM_TETRABORATE,
+    DEDX_MERCURIC_IODIDE,
+    DEDX_METHANE,
+    DEDX_METHANOL,
+    DEDX_MIX_D_WAX,
+    DEDX_MS20_TISSUE_SUBSTITUTE,
+    DEDX_MUSCLE_SKELETAL,
+    DEDX_MUSCLE_STRIATED,
+    DEDX_MUSCLE_EQUIVALENT_LIQUID_WITH_SUCROSE,
+    DEDX_MUSCLE_EQUIVALENT_LIQUID_WITHOUT_SUCROSE,
+    DEDX_NAPHTHALENE,
+    DEDX_NITROBENZENE,
+    DEDX_NITROUS_OXIDE,
+    DEDX_NYLON_DUPONT_ELVAMIDE_8062,
+    DEDX_NYLON_TYPE_6_AND_6_6,
+    DEDX_NYLON_TYPE_6_10,
+    DEDX_NYLON_TYPE_11_RILSAN,
+    DEDX_OCTANE_LIQUID,
+    DEDX_PARAFFIN_WAX,
+    DEDX_N_PENTANE,
+    DEDX_PHOTOGRAPHIC_EMULSION,
+    DEDX_PLASTIC_SCINTILLATOR_VINYLTOLUENE_BASED,
+    DEDX_PLUTONIUM_DIOXIDE,
+    DEDX_POLYACRYLONITRILE,
+    DEDX_POLYCARBONATE_MAKROLON_LEXAN,
+    DEDX_POLYCHLOROSTYRENE,
+    DEDX_POLYETHYLENE,
+    DEDX_MYLAR,
+    DEDX_LUCITE_PERSPEX_PMMA,
+    DEDX_POLYOXYMETHYLENE,
+    DEDX_POLYPROPYLENE,
+    DEDX_POLYSTYRENE,
+    DEDX_POLYTETRAFLUOROETHYLENE,
+    DEDX_POLYTRIFLUOROCHLOROETHYLENE,
+    DEDX_POLYVINYL_ACETATE,
+    DEDX_POLYVINYL_ALCOHOL,
+    DEDX_POLYVINYL_BUTYRAL,
+    DEDX_POLYVINYL_CHLORIDE,
+    DEDX_POLYVINYLIDENE_CHLORIDE_SARAN,
+    DEDX_POLYVINYLIDENE_FLUORIDE,
+    DEDX_POLYVINYL_PYRROLIDONE,
+    DEDX_POTASSIUM_IODIDE,
+    DEDX_POTASSIUM_OXIDE,
+    DEDX_PROPANE,
+    DEDX_PROPANE_LIQUID,
+    DEDX_N_PROPYL_ALCOHOL,
+    DEDX_PYRIDINE,
+    DEDX_RUBBER_BUTYL,
+    DEDX_RUBBER_NATURAL,
+    DEDX_RUBBER_NEOPRENE,
+    DEDX_SILICON_DIOXIDE,
+    DEDX_SILVER_BROMIDE,
+    DEDX_SILVER_CHLORIDE,
+    DEDX_SILVER_HALIDES_IN_PHOTOGRAPHIC_EMULSION,
+    DEDX_SILVER_IODIDE,
+    DEDX_SKIN_ICRP,
+    DEDX_SODIUM_CARBONATE,
+    DEDX_SODIUM_IODIDE,
+    DEDX_SODIUM_MONOXIDE,
+    DEDX_SODIUM_NITRATE,
+    DEDX_STILBENE,
+    DEDX_SUCROSE,
+    DEDX_TERPHENYL,
+    DEDX_TESTES_ICRP,
+    DEDX_TETRACHLOROETHYLENE,
+    DEDX_THALLIUM_CHLORIDE,
+    DEDX_TISSUE_SOFT_ICRP,
+    DEDX_TISSUE_SOFT_ICRU_FOUR_COMPONENT,
+    DEDX_TISSUE_EQUIVALENT_GAS_METHANE_BASED,
+    DEDX_TISSUE_EQUIVALENT_GAS_PROPANE_BASED,
+    DEDX_TITANIUM_DIOXIDE,
+    DEDX_TOLUENE,
+    DEDX_TRICHLOROETHYLENE,
+    DEDX_TRIETHYL_PHOSPHATE,
+    DEDX_TUNGSTEN_HEXAFLUORIDE,
+    DEDX_URANIUM_DICARBIDE,
+    DEDX_URANIUM_MONOCARBIDE,
+    DEDX_URANIUM_OXIDE,
+    DEDX_UREA,
+    DEDX_VALINE,
+    DEDX_VITON_FLUOROELASTOMER,
+    DEDX_WATER_LIQUID,
+    DEDX_WATER_VAPOR,
+    DEDX_XYLENE
+};
+
+/** @} */ /* ions_and_materials */
+
+/** @defgroup special_ions Special particle identifiers
+ *  @brief Identifiers for particles that do not map directly to atomic number.
+ *  @{
+ */
+#define DEDX_PROTON 1        /**< Proton (Z=1, same as DEDX_HYDROGEN) */
+#define DEDX_ELECTRON 1001   /**< Electron */
+#define DEDX_POSITRON 1002   /**< Positron */
+#define DEDX_PIMINUS 1003    /**< Pi minus */
+#define DEDX_PIPLUS 1004     /**< Pi plus */
+#define DEDX_PIZERO 1005     /**< Pi zero */
+#define DEDX_ANTIPROTON 1006 /**< Antiproton */
+/** @} */
+
+/** @defgroup aliases Common material aliases
+ *  @{
+ */
+#define DEDX_WATER DEDX_WATER_LIQUID
+#define DEDX_AIR DEDX_AIR_DRY_NEAR_SEA_LEVEL
+#define DEDX_PMMA DEDX_LUCITE_PERSPEX_PMMA
+#define DEDX_PERSPEX DEDX_LUCITE_PERSPEX_PMMA
+#define DEDX_LUCITE DEDX_LUCITE_PERSPEX_PMMA
+#define DEDX_TEFLON DEDX_POLYTETRAFLUOROETHYLENE
+#define DEDX_CONCRETE DEDX_CONCRETE_PORTLAND
+#define DEDX_CAESIUM DEDX_CESIUM
+/** @} */
+
+#endif // DEDX_ELEMENTS_H
diff --git a/libdedx/dedx_error.h b/libdedx/dedx_error.h
index a8784e4..a7d2d81 100644
--- a/libdedx/dedx_error.h
+++ b/libdedx/dedx_error.h
@@ -1,19 +1,19 @@
-#ifndef DEDX_ERROR_H_INCLUDED
-#define DEDX_ERROR_H_INCLUDED
+#ifndef DEDX_ERROR_H
+#define DEDX_ERROR_H
 
 /**
  * @file dedx_error.h
  * @brief Error codes returned via the int *err output parameter.
  *
- * Functions signal success by setting *err = DEDX_OK (0).
- * On failure they set *err to one of the values below and, unless
+ * Functions signal success by setting `*err = DEDX_OK`.
+ * On failure they set `*err` to one of the values below and, unless
  * stated otherwise, return 0 / NULL / -1.
  */
 
 /** No error. */
 #define DEDX_OK 0
 
-/** @defgroup err_file File I/O errors (1–11)
+/** @defgroup err_file File I/O errors (1–10)
  * @{ */
 #define DEDX_ERR_NO_COMPOS_FILE 1   /**< compos.txt not found */
 #define DEDX_ERR_NO_GAS_FILE 2      /**< gas_states.dat not found */
@@ -23,9 +23,8 @@
 #define DEDX_ERR_WRITE_FAILED 6     /**< cannot write to disk */
 #define DEDX_ERR_NO_ENERGY_FILE 7   /**< energy .dat source file unreadable */
 #define DEDX_ERR_NO_DATA_FILE 8     /**< stopping-power .dat source file unreadable */
-#define DEDX_ERR_NO_NAMES_FILE 9    /**< short_names file unreadable (reserved) */
+#define DEDX_ERR_NO_NAMES_FILE 9    /**< short_names file unreadable (reserved internal code) */
 #define DEDX_ERR_NO_COMPOSITION 10  /**< elemental composition file unreadable */
-#define DEDX_ERR_NO_ATIMA_FILE 11   /**< ATIMA composition file unreadable (reserved) */
 /** @} */
 
 /** @defgroup err_bounds Bounds errors (101)
@@ -40,11 +39,12 @@
 #define DEDX_ERR_INVALID_DATASET_ID 203      /**< dataset ID does not exist in workspace */
 #define DEDX_ERR_NOT_AN_ELEMENT 204          /**< ID does not correspond to an atomic element */
 #define DEDX_ERR_ESTAR_NOT_IMPL 205          /**< ESTAR program is not implemented */
-#define DEDX_ERR_ION_NOT_SUPPORTED_MSTAR 206 /**< ion not supported for MSTAR (reserved) */
+#define DEDX_ERR_ION_NOT_SUPPORTED_MSTAR 206 /**< ion not supported for MSTAR (reserved legacy code) */
 #define DEDX_ERR_ION_NOT_SUPPORTED 207       /**< ion not supported by requested program */
 #define DEDX_ERR_RHO_REQUIRED 208            /**< target density rho must be provided */
 #define DEDX_ERR_ION_A_REQUIRED 209          /**< nucleon number ion_a must be provided */
 #define DEDX_ERR_INVALID_I_VALUE 210         /**< mean excitation potential must be > 0 */
+#define DEDX_ERR_INCONSISTENT_COMPOUND 211   /**< inconsistent compound specification */
 /** @} */
 
 /** @defgroup err_memory Memory errors (301)
@@ -52,4 +52,4 @@
 #define DEDX_ERR_NO_MEMORY 301 /**< memory allocation failed */
 /** @} */
 
-#endif /* DEDX_ERROR_H_INCLUDED */
+#endif /* DEDX_ERROR_H */
diff --git a/libdedx/dedx_file.c b/libdedx/dedx_file.c
index e5dfdd3..d56ac77 100644
--- a/libdedx/dedx_file.c
+++ b/libdedx/dedx_file.c
@@ -16,71 +16,50 @@
 */
 #include "dedx_file.h"
 
-char *_dedx_get_program_file(int program) {
-    char *path;
+const char *dedx_internal_get_program_file(int program) {
     switch (program) {
-    case 1:
-        path = "ASTAR";
-        break;
-    case 2:
-        path = "PSTAR";
-        break;
-    case 3:
-        path = "ESTAR";
-        break;
-    case 4:
-        path = "MSTAR";
-        break;
-    case 5:
-        path = "ICRU73";
-        break;
-    case 6:
-        path = "ICRU73_NEW";
-        break;
-    case 7:
-        path = "ICRU_ASTAR";
-        break;
-    case 8:
-        path = "ICRU_PSTAR";
-        break;
+    case DEDX_ASTAR:
+        return "ASTAR";
+    case DEDX_PSTAR:
+        return "PSTAR";
+    case DEDX_ESTAR:
+        return "ESTAR";
+    case DEDX_MSTAR:
+        return "MSTAR";
+    case DEDX_ICRU73_OLD:
+        return "ICRU73";
+    case DEDX_ICRU73:
+        return "ICRU73_NEW";
+    case DEDX_ICRU49:
+        return "ICRU_ASTAR";
+    case _DEDX_0008:
+        return "ICRU_PSTAR";
     default:
-        path = "";
+        return "";
     }
-    return path;
 }
 
-char *_dedx_get_energy_file(int program) {
-    char *path;
+const char *dedx_internal_get_energy_file(int program) {
     switch (program) {
-    case 1:
-        path = "astarEng";
-        break;
-    case 2:
-        path = "pstarEng";
-        break;
-    case 3:
-        path = "estarEng";
-        break;
-    case 4:
-        path = "mstarEng";
-        break;
-    case 5:
-        path = "icru73Eng";
-        break;
-    case 6:
-        path = "icru73_newEng";
-        break;
-    case 7:
-        path = "icru_astarEng";
-        break;
-    case 8:
-        path = "icru_pstarEng";
-        break;
-    case 101:
-        path = "betheEng";
-        break;
+    case DEDX_ASTAR:
+        return "astarEng";
+    case DEDX_PSTAR:
+        return "pstarEng";
+    case DEDX_ESTAR:
+        return "estarEng";
+    case DEDX_MSTAR:
+        return "mstarEng";
+    case DEDX_ICRU73_OLD:
+        return "icru73Eng";
+    case DEDX_ICRU73:
+        return "icru73_newEng";
+    case DEDX_ICRU49:
+        return "icru_astarEng";
+    case _DEDX_0008:
+        return "icru_pstarEng";
+    case DEDX_BETHE_EXT00:
+        return "betheEng";
     default:
-        path = "";
+        return "";
     }
-    return path;
 }
diff --git a/libdedx/dedx_file.h b/libdedx/dedx_file.h
index d14220b..c9bf4f3 100644
--- a/libdedx/dedx_file.h
+++ b/libdedx/dedx_file.h
@@ -1,8 +1,10 @@
-#ifndef DEDX_FILE_H_INCLUDED
-#define DEDX_FILE_H_INCLUDED
+#ifndef DEDX_FILE_H
+#define DEDX_FILE_H
 
-char *_dedx_get_program_file(int program);
+#include "dedx.h"
 
-char *_dedx_get_energy_file(int program);
+const char *dedx_internal_get_program_file(int program);
 
-#endif // DEDX_FILE_H_INCLUDED
+const char *dedx_internal_get_energy_file(int program);
+
+#endif // DEDX_FILE_H
diff --git a/libdedx/dedx_file_access.c b/libdedx/dedx_file_access.c
index 022aef3..0941bfa 100644
--- a/libdedx/dedx_file_access.c
+++ b/libdedx/dedx_file_access.c
@@ -17,6 +17,17 @@
 
 #include "dedx_file_access.h"
 
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <sys/stat.h>
+
+#include "dedx_config.h"
+#include "dedx_file.h"
+#include "dedx_split.h"
+
+#define DEDX_PATH_SIZE 512
+
 /* Resolve the data directory once per process.
  * Prefers DEDX_DATA_PATH_LOCAL (build tree) if it contains the sentinel file,
  * otherwise falls back to DEDX_DATA_PATH (install prefix).
@@ -26,7 +37,7 @@
  * concurrent first calls. The library as a whole is not thread-safe; this
  * function should be fixed with call_once() as part of a broader thread-safety
  * effort. See https://github.com/APTG/libdedx/issues (thread-safety issue). */
-static const char *_dedx_get_data_path(void) {
+static const char *get_data_path(void) {
     static char resolved[DEDX_PATH_SIZE];
     static int done = 0;
     struct stat st;
@@ -34,7 +45,7 @@ static const char *_dedx_get_data_path(void) {
 
     if (!done) {
         snprintf(probe, sizeof(probe), "%s%s", DEDX_DATA_PATH_LOCAL, "icru_pstarEng.dat");
-        if (stat(probe, &st) == 0)
+        if (DEDX_DATA_PATH_LOCAL[0] != '\0' && stat(probe, &st) == 0)
             snprintf(resolved, sizeof(resolved), "%s", DEDX_DATA_PATH_LOCAL);
         else
             snprintf(resolved, sizeof(resolved), "%s", DEDX_DATA_PATH); /* LCOV_EXCL_LINE */
@@ -43,14 +54,14 @@ static const char *_dedx_get_data_path(void) {
     return resolved;
 }
 
-void _dedx_convert_to_binary(char *path, char *output, int *err) {
+static void convert_to_binary(char *path, char *output, int *err) {
     FILE *fp;
     FILE *out;
     char line[100];
     int datalines = 0;
     int i;
     unsigned int length;
-    float data[_DEDX_MAXELEMENTS];
+    float data[DEDX_MAX_ELEMENTS];
     char **temp = NULL;
     stopping_data container;
 
@@ -68,34 +79,44 @@ void _dedx_convert_to_binary(char *path, char *output, int *err) {
         }
         return;
     } /* LCOV_EXCL_STOP */
-    while (!feof(fp)) {
-        if (fgets(line, 100, fp) == NULL) {
-        }
+    while (fgets(line, sizeof(line), fp) != NULL) {
         if (line[0] == '#') {
             length = 0;
             i = 0;
-            memset(&data, 0, _DEDX_MAXELEMENTS);
-            temp = _dedx_split(line, ':', &length, 100);
+            memset(data, 0, sizeof(data));
+            temp = dedx_internal_split(line, ':', &length, 100);
+            if (temp == NULL) {
+                *err = DEDX_ERR_NO_MEMORY;
+                break;
+            }
             temp[0][0] = ' ';
             datalines = atoi(temp[2]);
-            while (i++ < datalines && !feof(fp)) {
-                if (fgets(line, 100, fp) == NULL) {
+            while (i++ < datalines) {
+                if (fgets(line, sizeof(line), fp) == NULL) {
+                    *err = DEDX_ERR_NO_DATA_FILE;
+                    dedx_internal_free_split_temp(temp, 10);
+                    temp = NULL;
+                    break;
                 }
                 data[i - 1] = atof(line);
             }
+            if (*err != DEDX_OK)
+                break;
             container.target = atoi(temp[0]);
             container.ion = atoi(temp[1]);
             container.length = datalines;
-            memcpy(&container.data, &data, _DEDX_MAXELEMENTS * sizeof(float));
+            memcpy(container.data, data, sizeof(data));
             fwrite(&container, sizeof(container), 1, out);
+            dedx_internal_free_split_temp(temp, 10);
+            temp = NULL;
         }
     }
-    free(temp);
+    dedx_internal_free_split_temp(temp, 10);
     fclose(fp);
     fclose(out);
 }
 
-void _dedx_read_binary_data(stopping_data *data, int prog, int ion, int target, int *err) {
+void dedx_internal_read_binary_data(stopping_data *data, int prog, int ion, int target, int *err) {
     const char *folder;
     char path[DEDX_PATH_SIZE];
     char input_path[DEDX_PATH_SIZE];
@@ -103,13 +124,13 @@ void _dedx_read_binary_data(stopping_data *data, int prog, int ion, int target,
     stopping_data dat;
 
     *err = DEDX_OK;
-    folder = _dedx_get_data_path();
-    snprintf(path, sizeof(path), "%s%s.bin", folder, _dedx_get_program_file(prog));
+    folder = get_data_path();
+    snprintf(path, sizeof(path), "%s%s.bin", folder, dedx_internal_get_program_file(prog));
 
     fp = fopen(path, "rb");
     if (fp == NULL) {
-        snprintf(input_path, sizeof(input_path), "%s%s.dat", folder, _dedx_get_program_file(prog));
-        _dedx_convert_to_binary(input_path, path, err);
+        snprintf(input_path, sizeof(input_path), "%s%s.dat", folder, dedx_internal_get_program_file(prog));
+        convert_to_binary(input_path, path, err);
         if (*err != DEDX_OK) {
             *err = DEDX_ERR_NO_BINARY_DATA;
             return;
@@ -130,17 +151,16 @@ void _dedx_read_binary_data(stopping_data *data, int prog, int ion, int target,
     fclose(fp);
 }
 
-void _dedx_convert_energy_binary(char *path, char *output, int *err) {
+static void convert_energy_binary(char *path, char *output, int *err) {
     FILE *fp;
     FILE *out;
     char line[100];
     int datalines;
     int i;
-    float data[_DEDX_MAXELEMENTS];
+    float data[DEDX_MAX_ELEMENTS];
 
     *err = DEDX_OK;
     fp = fopen(path, "r");
-    // TODO: Next line wont work after installation, unless user is root.
     out = fopen(output, "wb+");
     if (fp == NULL || out == NULL) { /* LCOV_EXCL_START */
         if (out == NULL) {
@@ -154,13 +174,23 @@ void _dedx_convert_energy_binary(char *path, char *output, int *err) {
         return;
     } /* LCOV_EXCL_STOP */
 
-    if (fgets(line, 100, fp) == NULL) {
+    memset(data, 0, sizeof(data));
+
+    if (fgets(line, sizeof(line), fp) == NULL) {
+        *err = DEDX_ERR_NO_ENERGY_FILE;
+        fclose(fp);
+        fclose(out);
+        return;
     }
     datalines = atoi(line);
 
     for (i = 0; i < datalines; i++) {
-        if (fgets(line, 100, fp) != NULL) {
-        };
+        if (fgets(line, sizeof(line), fp) == NULL) {
+            *err = DEDX_ERR_NO_ENERGY_FILE;
+            fclose(fp);
+            fclose(out);
+            return;
+        }
         data[i] = atof(line);
     }
     fwrite(&data, sizeof(data), 1, out);
@@ -168,20 +198,20 @@ void _dedx_convert_energy_binary(char *path, char *output, int *err) {
     fclose(out);
 }
 
-void _dedx_read_energy_data(float *energy, int prog, int *err) {
+void dedx_internal_read_energy_data(float *energy, int prog, int *err) {
     const char *folder;
     char path[DEDX_PATH_SIZE];
     char input_path[DEDX_PATH_SIZE];
     FILE *fp;
 
     *err = DEDX_OK;
-    folder = _dedx_get_data_path();
-    snprintf(path, sizeof(path), "%s%s.bin", folder, _dedx_get_energy_file(prog));
+    folder = get_data_path();
+    snprintf(path, sizeof(path), "%s%s.bin", folder, dedx_internal_get_energy_file(prog));
 
     fp = fopen(path, "rb");
     if (fp == NULL) {
-        snprintf(input_path, sizeof(input_path), "%s%s.dat", folder, _dedx_get_energy_file(prog));
-        _dedx_convert_energy_binary(input_path, path, err);
+        snprintf(input_path, sizeof(input_path), "%s%s.dat", folder, dedx_internal_get_energy_file(prog));
+        convert_energy_binary(input_path, path, err);
         if (*err != DEDX_OK) {
             *err = DEDX_ERR_NO_BINARY_ENERGY;
             return;
@@ -192,12 +222,12 @@ void _dedx_read_energy_data(float *energy, int prog, int *err) {
             return;
         }
     }
-    if (fread(energy, sizeof(float) * _DEDX_MAXELEMENTS, 1, fp) == 0) {
+    if (fread(energy, sizeof(float) * DEDX_MAX_ELEMENTS, 1, fp) == 0) {
     }
     fclose(fp);
 }
 
-float _dedx_read_effective_charge(int id, int *err) {
+float dedx_internal_read_effective_charge(int id, int *err) {
     const char *folder;
     char line[100];
     char **temp;
@@ -210,7 +240,7 @@ float _dedx_read_effective_charge(int id, int *err) {
     if (id < 99) {
         return (float) id;
     }
-    folder = _dedx_get_data_path();
+    folder = get_data_path();
     items = 0;
     temp = NULL;
     snprintf(path, sizeof(path), "%s%s", folder, "effective_charge.dat");
@@ -221,22 +251,24 @@ float _dedx_read_effective_charge(int id, int *err) {
         return 0;
     }
 
-    while (!feof(fp)) {
-        if (fgets(line, 100, fp) != NULL) {
+    while (fgets(line, sizeof(line), fp) != NULL) {
+        temp = dedx_internal_split(line, '\t', &items, 100);
+        if (temp == NULL) {
+            *err = DEDX_ERR_NO_MEMORY;
+            break;
         }
-        temp = _dedx_split(line, '\t', &items, 100);
         if (atoi(temp[0]) == id) {
             charge = atof(temp[1]);
-            _dedx_free_split_temp(temp);
+            dedx_internal_free_split_temp(temp, 10);
             break;
         }
-        _dedx_free_split_temp(temp);
+        dedx_internal_free_split_temp(temp, 10);
     }
     fclose(fp);
     return charge;
 }
 
-size_t _dedx_target_is_gas(int target, int *err) {
+size_t dedx_internal_target_is_gas(int target, int *err) {
     const char *folder;
     size_t is_gas;
     char str[100];
@@ -245,7 +277,7 @@ size_t _dedx_target_is_gas(int target, int *err) {
 
     *err = DEDX_OK;
     is_gas = 0;
-    folder = _dedx_get_data_path();
+    folder = get_data_path();
     snprintf(path, sizeof(path), "%s%s", folder, "gas_states.dat");
     fp = fopen(path, "r");
     if (fp == NULL) {
@@ -253,9 +285,7 @@ size_t _dedx_target_is_gas(int target, int *err) {
         return 0;
     }
 
-    while (!feof(fp)) {
-        if (fgets(str, 100, fp) == NULL) {
-        }
+    while (fgets(str, sizeof(str), fp) != NULL) {
         if (atoi(str) == target)
             is_gas = 1;
     }
@@ -263,7 +293,7 @@ size_t _dedx_target_is_gas(int target, int *err) {
     return is_gas;
 }
 
-float _dedx_read_density(int id, int *err) {
+float dedx_internal_read_density(int id, int *err) {
     const char *folder;
     float density;
     char str[100];
@@ -275,7 +305,7 @@ float _dedx_read_density(int id, int *err) {
     *err = DEDX_OK;
     density = 1.0;
     items = 0;
-    folder = _dedx_get_data_path();
+    folder = get_data_path();
     snprintf(path, sizeof(path), "%s%s", folder, "compos.txt");
     fp = fopen(path, "r");
 
@@ -284,16 +314,18 @@ float _dedx_read_density(int id, int *err) {
         return 0;
     }
 
-    while (!feof(fp)) {
-        if (fgets(str, 100, fp) == NULL) {
+    while (fgets(str, sizeof(str), fp) != NULL) {
+        temp = dedx_internal_split(str, '\t', &items, 100);
+        if (temp == NULL) {
+            *err = DEDX_ERR_NO_MEMORY;
+            break;
         }
-        temp = _dedx_split(str, '\t', &items, 100);
         if (atoi(temp[0]) == id) {
             density = atof(temp[1]);
-            _dedx_free_split_temp(temp);
+            dedx_internal_free_split_temp(temp, 10);
             break;
         }
-        _dedx_free_split_temp(temp);
+        dedx_internal_free_split_temp(temp, 10);
     }
     fclose(fp);
     if (density == 0.0) {
@@ -302,7 +334,7 @@ float _dedx_read_density(int id, int *err) {
     return density;
 }
 
-float _dedx_get_i_value(int target, int state, int *err) {
+float dedx_internal_get_i_value(int target, int state, int *err) {
     const char *folder;
     float pot;
     char path[DEDX_PATH_SIZE];
@@ -311,9 +343,10 @@ float _dedx_get_i_value(int target, int state, int *err) {
     unsigned int items;
     FILE *fp;
 
+    *err = DEDX_OK;
     pot = 0.0;
     items = 0;
-    folder = _dedx_get_data_path();
+    folder = get_data_path();
     snprintf(path, sizeof(path), "%s%s", folder, "compos.txt");
 
     fp = fopen(path, "r");
@@ -321,10 +354,12 @@ float _dedx_get_i_value(int target, int state, int *err) {
         *err = DEDX_ERR_NO_COMPOS_FILE;
         return 0;
     }
-    while (!feof(fp)) {
-        if (fgets(str, 100, fp) == NULL) {
+    while (fgets(str, sizeof(str), fp) != NULL) {
+        temp = dedx_internal_split(str, '\t', &items, 100);
+        if (temp == NULL) {
+            *err = DEDX_ERR_NO_MEMORY;
+            break;
         }
-        temp = _dedx_split(str, '\t', &items, 100);
         if (atoi(temp[0]) == target) {
             if (items == 4) {
                 if (state == atoi(temp[3])) {
@@ -332,11 +367,15 @@ float _dedx_get_i_value(int target, int state, int *err) {
                 }
             } else {
                 pot = atof(temp[2]);
-                if (state == 2 && _dedx_target_is_gas(target, err) == 0 && target <= 99)
+                if (state == 2 && dedx_internal_target_is_gas(target, err) == 0 && target <= 99)
                     pot = pot * 1.13;
+                if (*err != DEDX_OK) {
+                    dedx_internal_free_split_temp(temp, 10);
+                    break;
+                }
             }
         }
-        _dedx_free_split_temp(temp);
+        dedx_internal_free_split_temp(temp, 10);
     }
     fclose(fp);
     if (pot == 0.0) {
@@ -345,7 +384,7 @@ float _dedx_get_i_value(int target, int state, int *err) {
     return pot;
 }
 
-void _dedx_get_composition(int target, float composition[][2], unsigned int *length, int *err) {
+void dedx_internal_get_composition(int target, float composition[][2], unsigned int *length, int *err) {
     const char *folder;
     char path[DEDX_PATH_SIZE];
     char str[100];
@@ -356,7 +395,7 @@ void _dedx_get_composition(int target, float composition[][2], unsigned int *len
 
     *err = DEDX_OK;
     *length = 0;
-    folder = _dedx_get_data_path();
+    folder = get_data_path();
     snprintf(path, sizeof(path), "%s%s", folder, "composition");
 
     fp = fopen(path, "r");
@@ -364,59 +403,26 @@ void _dedx_get_composition(int target, float composition[][2], unsigned int *len
         *err = DEDX_ERR_NO_COMPOSITION;
         return;
     }
-    while (!feof(fp)) {
-        if (fgets(str, 100, fp) != NULL)
-            if (str[0] == '#') {
-                if (atoi(&str[1]) == target) {
-                    while (!feof(fp)) {
-                        if (fgets(str, 100, fp) != NULL) {
-                        }
-                        f_temp = atof(str);
-                        if (f_temp == 0.0)
-                            break;
-                        temp = _dedx_split(str, ':', &items, 100);
-                        composition[(*length)][0] = atof(temp[0]);
-                        composition[(*length)++][1] = atof(temp[1]);
-                        _dedx_free_split_temp(temp);
-                    }
-
-                    break;
-                }
+    while (fgets(str, sizeof(str), fp) != NULL) {
+        if (str[0] != '#')
+            continue;
+        if (atoi(&str[1]) != target)
+            continue;
+
+        while (fgets(str, sizeof(str), fp) != NULL) {
+            f_temp = atof(str);
+            if (f_temp == 0.0)
+                break;
+            temp = dedx_internal_split(str, ':', &items, 100);
+            if (temp == NULL) {
+                *err = DEDX_ERR_NO_MEMORY;
+                break;
             }
-    }
-    fclose(fp);
-}
-
-float *_dedx_get_atima_data(int target, int *err) {
-    const char *folder;
-    char path[DEDX_PATH_SIZE];
-    char str[100];
-    float *compos;
-    unsigned int items;
-    char **temp;
-    FILE *fp;
-
-    *err = DEDX_OK;
-    folder = _dedx_get_data_path();
-    snprintf(path, sizeof(path), "%s%s", folder, "atima_compos");
-
-    fp = fopen(path, "r");
-    if (fp == NULL) {
-        *err = DEDX_ERR_NO_COMPOSITION;
-        return NULL;
-    }
-    compos = (float *) malloc(sizeof(float) * 4);
-    while (!feof(fp) && fgets(str, 100, fp) != NULL) {
-        temp = _dedx_split(str, '\t', &items, 100);
-        if (atoi(temp[0]) == target) {
-            compos = (float *) malloc(sizeof(float) * 4);
-            compos[0] = atof(temp[1]);
-            compos[1] = atof(temp[2]);
-            compos[2] = atof(temp[3]);
-            compos[3] = atof(temp[4]);
+            composition[(*length)][0] = atof(temp[0]);
+            composition[(*length)++][1] = atof(temp[1]);
+            dedx_internal_free_split_temp(temp, 10);
         }
-        _dedx_free_split_temp(temp);
+        break;
     }
     fclose(fp);
-    return compos;
 }
diff --git a/libdedx/dedx_file_access.h b/libdedx/dedx_file_access.h
index b51deac..775c180 100644
--- a/libdedx/dedx_file_access.h
+++ b/libdedx/dedx_file_access.h
@@ -1,32 +1,60 @@
-#ifndef DEDX_FILE_ACCESS_H_INCLUDED
-#define DEDX_FILE_ACCESS_H_INCLUDED
-#include <stdio.h>
-#include <stdlib.h>
-#include <string.h>
-
-#include "dedx.h"
-// #include "dedx_constants.h"
-#include "dedx_file.h"
+#ifndef DEDX_FILE_ACCESS_H
+#define DEDX_FILE_ACCESS_H
+
+#include <stddef.h>
+
 #include "dedx_stopping_data.h"
-// #include "utils/dedx_string.h"
-#include <sys/stat.h>
 
-#include "dedx_config.h"
-#include "dedx_split.h"
+/** @brief Read one tabulated stopping-power dataset from the cached binary file.
+ *  @param[out] data    Destination dataset buffer.
+ *  @param[in]  prog    Program identifier.
+ *  @param[in]  ion     Projectile identifier.
+ *  @param[in]  target  Target identifier.
+ *  @param[out] err     Set to DEDX_OK or a file/data error code.
+ */
+void dedx_internal_read_binary_data(stopping_data *data, int prog, int ion, int target, int *err);
+
+/** @brief Read the energy grid associated with a tabulated program.
+ *  @param[out] energy  Destination array of DEDX_MAX_ELEMENTS values.
+ *  @param[in]  prog    Program identifier.
+ *  @param[out] err     Set to DEDX_OK or a file/data error code.
+ */
+void dedx_internal_read_energy_data(float *energy, int prog, int *err);
+
+/** @brief Read an effective charge override for a heavy target or ion ID.
+ *  @param[in]  id   Element or material identifier.
+ *  @param[out] err  Set to DEDX_OK or a file/data error code.
+ *  @return Effective charge, or @p id itself for elemental IDs below 99.
+ */
+float dedx_internal_read_effective_charge(int id, int *err);
 
-#define DEDX_PATH_SIZE 512
+/** @brief Check whether a target is listed as gaseous in gas_states.dat.
+ *  @param[in]  target  Target identifier.
+ *  @param[out] err     Set to DEDX_OK or a file/data error code.
+ *  @return Non-zero if the target is gaseous, zero otherwise.
+ */
+size_t dedx_internal_target_is_gas(int target, int *err);
 
-void _dedx_convert_to_binary(char *path, char *output, int *err);
-void _dedx_read_binary_data(stopping_data *data, int prog, int ion, int target, int *err);
+/** @brief Read the default density for a target from compos.txt.
+ *  @param[in]  id   Target identifier.
+ *  @param[out] err  Set to DEDX_OK or a file/data error code.
+ *  @return Density in g/cm^3, or 0 on failure.
+ */
+float dedx_internal_read_density(int id, int *err);
 
-void _dedx_convert_energy_binary(char *path, char *output, int *err);
-void _dedx_read_energy_data(float *energy, int prog, int *err);
-float _dedx_read_effective_charge(int id, int *err);
-size_t _dedx_target_is_gas(int target, int *err);
-float _dedx_read_density(int id, int *err);
-// void _dedx_get_short_name(int id,char * name, int *err);
-float _dedx_get_i_value(int target, int state, int *err);
-void _dedx_get_composition(int target, float composition[][2], unsigned int *length, int *err);
+/** @brief Read the mean excitation potential for a target and state.
+ *  @param[in]  target  Target identifier.
+ *  @param[in]  state   Compound-state selector used by compos.txt.
+ *  @param[out] err     Set to DEDX_OK or a file/data error code.
+ *  @return Mean excitation potential in eV, or 0 on failure.
+ */
+float dedx_internal_get_i_value(int target, int state, int *err);
 
-float *_dedx_get_atima_data(int target, int *err);
-#endif // DEDX_FILE_ACCESS_H_INCLUDED
+/** @brief Read the elemental composition of a compound target.
+ *  @param[in]  target       Compound target identifier.
+ *  @param[out] composition  Output array of `(element_id, weight)` pairs.
+ *  @param[out] length       Number of populated rows in @p composition.
+ *  @param[out] err          Set to DEDX_OK or a file/data error code.
+ */
+void dedx_internal_get_composition(int target, float composition[][2], unsigned int *length, int *err);
+#endif // DEDX_FILE_ACCESS_H
diff --git a/libdedx/dedx_gold_struct.h b/libdedx/dedx_gold_struct.h
deleted file mode 100644
index 23fa846..0000000
--- a/libdedx/dedx_gold_struct.h
+++ /dev/null
@@ -1,17 +0,0 @@
-#ifndef DEDX_GOLD_STRUCT_H_INCLUDED
-#define DEDX_GOLD_STRUCT_H_INCLUDED
-
-/*typedef struct
-{
-    double e_min;
-    double e_max;
-    double epsilon;
-    double h;
-    double e_zero;
-    double e_extr;
-    double f_extr;
-    double e_sewn;
-    double f_sewn;
-} gold_struct;*/
-
-#endif // DEDX_GOLD_STRUCT_H_INCLUDED
diff --git a/libdedx/dedx_lookup_accelerator.h b/libdedx/dedx_lookup_accelerator.h
index b6a8216..ebb3337 100644
--- a/libdedx/dedx_lookup_accelerator.h
+++ b/libdedx/dedx_lookup_accelerator.h
@@ -1,13 +1,16 @@
-#ifndef DEDX_LOOKUP_ACCELERATOR_H_INCLUDED
-#define DEDX_LOOKUP_ACCELERATOR_H_INCLUDED
+#ifndef DEDX_LOOKUP_ACCELERATOR_H
+#define DEDX_LOOKUP_ACCELERATOR_H
 
-// #include "dedx.h"
+/** @brief Cached interval index for repeated spline evaluations.
+ *
+ *  The spline evaluator stores the last successful interval in @c cache so
+ *  nearby lookups can skip the binary search. The hit/miss counters are
+ *  currently informational only.
+ */
+typedef struct {
+    int cache;
+    int hits;
+    int miss;
+} _dedx_lookup_accelerator;
 
-/* typedef struct */
-/* { */
-/*     int cache; */
-/*     int hits; */
-/*     int miss; */
-/* } lookup_accelerator; */
-
-#endif // DEDX_LOOKUP_ACCELERATOR_H_INCLUDED
+#endif // DEDX_LOOKUP_ACCELERATOR_H
diff --git a/libdedx/dedx_lookup_data.h b/libdedx/dedx_lookup_data.h
index 95c23ba..4eb88d8 100644
--- a/libdedx/dedx_lookup_data.h
+++ b/libdedx/dedx_lookup_data.h
@@ -1,17 +1,23 @@
-#ifndef DEDX_LOOKUP_DATA_H_INCLUDED
-#define DEDX_LOOKUP_DATA_H_INCLUDED
+#ifndef DEDX_LOOKUP_DATA_H
+#define DEDX_LOOKUP_DATA_H
 
-// #include "dedx.h"
+#include "dedx_elements.h"
+#include "dedx_lookup_accelerator.h"
+#include "dedx_spline_base.h"
 
-/* typedef struct */
-/* { */
-/*     spline_base base[MAXELEMENTS]; */
-/*     int n; */
-/*     int prog; */
-/*     int target; */
-/*     int ion; */
-/*     int datapoints; */
-/*     lookup_accelerator acc; */
-/* } lookup_data; */
+/** @brief Internal representation of one loaded stopping-power table.
+ *
+ *  Each dataset stores the spline coefficients, source program id, ion/target
+ *  identifiers, and the lookup accelerator used by repeated spline queries.
+ */
+typedef struct _dedx_lookup_data {
+    _dedx_spline_base base[DEDX_MAX_ELEMENTS];
+    int n;
+    int prog;
+    int target;
+    int ion;
+    int datapoints;
+    _dedx_lookup_accelerator acc;
+} _dedx_lookup_data;
 
-#endif // DEDX_LOOKUP_DATA_H_INCLUDED
+#endif // DEDX_LOOKUP_DATA_H
diff --git a/libdedx/dedx_mpaul.c b/libdedx/dedx_mpaul.c
index a2858ae..c11e5fd 100644
--- a/libdedx/dedx_mpaul.c
+++ b/libdedx/dedx_mpaul.c
@@ -16,13 +16,13 @@
 */
 #include "dedx_mpaul.h"
 
-float _dedx_calculate_mspaul_coef(char mode, int ion, int target, float energy) {
+float dedx_internal_calculate_mspaul_coef(char mode, int ion, int target, float energy) {
     if (ion == 2)
         return 1;
     float output;
     int err = 0;
     float f;
-    f = _dedx_read_effective_charge(target, &err);
+    f = dedx_internal_read_effective_charge(target, &err);
 
     double a = 5.0;
     double b = -1;
diff --git a/libdedx/dedx_mpaul.h b/libdedx/dedx_mpaul.h
index e9f8149..ef1f730 100644
--- a/libdedx/dedx_mpaul.h
+++ b/libdedx/dedx_mpaul.h
@@ -1,9 +1,9 @@
-#ifndef DEDX_MPAUL_H_INCLUDED
-#define DEDX_MPAUL_H_INCLUDED
+#ifndef DEDX_MPAUL_H
+#define DEDX_MPAUL_H
 #include <math.h>
 
 #include "dedx_file_access.h"
 
-float _dedx_calculate_mspaul_coef(char mode, int ion, int target, float energy);
+float dedx_internal_calculate_mspaul_coef(char mode, int ion, int target, float energy);
 
-#endif // DEDX_MPAUL_H_INCLUDED
+#endif // DEDX_MPAUL_H
diff --git a/libdedx/dedx_mstar.c b/libdedx/dedx_mstar.c
index 9bbf2dc..7fdc4c5 100644
--- a/libdedx/dedx_mstar.c
+++ b/libdedx/dedx_mstar.c
@@ -16,46 +16,48 @@
 */
 #include "dedx_mstar.h"
 
-void _evaluate_compound_state_mstar(dedx_config *config, int *err) { /* LCOV_EXCL_START */
-    if (config->compound_state == DEDX_GAS) {
-        if (config->mstar_mode == 'a')
-            config->mstar_mode = 'g';
-        else
-            config->mstar_mode = 'h';
-    } else if (config->compound_state == DEDX_CONDENSED) {
-        if (config->mstar_mode == 'a')
-            config->mstar_mode = 'c';
-        else
-            config->mstar_mode = 'd';
-    }
-    *err = DEDX_OK;
-} /* LCOV_EXCL_STOP */
+#include "dedx_file_access.h"
 
-void _dedx_convert_energy_to_mstar(
-    stopping_data *in, stopping_data *out, char state, dedx_config *config, float *energy) {
-    //	int err = 0;
-    if (state == 'a' || state == 'b') {
-        if (config->compound_state) {
-            if (state == 'a')
-                state = 'g';
-            else
-                state = 'h';
+static char resolve_mstar_mode(char state, dedx_config *config, int *err) {
+    int target_state = config->compound_state;
 
+    if (state != DEDX_MSTAR_MODE_A && state != DEDX_MSTAR_MODE_B)
+        return state;
+
+    if (target_state == DEDX_DEFAULT_STATE) {
+        if (dedx_internal_target_is_gas(config->target, err) != 0) {
+            target_state = DEDX_GAS;
         } else {
-            if (state == 'a') {
-                state = 'c';
-            } else {
-                state = 'd';
-            }
+            if (*err != DEDX_OK)
+                return state;
+            target_state = DEDX_CONDENSED;
         }
     }
+
+    if (target_state == DEDX_GAS) {
+        return (state == DEDX_MSTAR_MODE_A) ? DEDX_MSTAR_MODE_G : DEDX_MSTAR_MODE_H;
+    }
+    return (state == DEDX_MSTAR_MODE_A) ? DEDX_MSTAR_MODE_C : DEDX_MSTAR_MODE_D;
+}
+
+void dedx_internal_evaluate_compound_state_mstar(dedx_config *config, int *err) { /* LCOV_EXCL_START */
+    *err = DEDX_OK;
+    config->mstar_mode = resolve_mstar_mode(config->mstar_mode, config, err);
+} /* LCOV_EXCL_STOP */
+
+void dedx_internal_convert_energy_to_mstar(
+    stopping_data *in, stopping_data *out, char state, dedx_config *config, float *energy, int *err) {
+    *err = DEDX_OK;
+    state = resolve_mstar_mode(state, config, err);
+    if (*err != DEDX_OK)
+        return;
     int n = in->length;
     int i;
     for (i = 0; i < n; i++) {
         energy[i] = energy[i] / 4.0;
     }
     for (i = 0; i < n; i++) {
-        out->data[i] = _dedx_calculate_mspaul_coef(state, in->ion, in->target, energy[i]) * in->data[i] * 1000;
+        out->data[i] = dedx_internal_calculate_mspaul_coef(state, in->ion, in->target, energy[i]) * in->data[i] * 1000;
     }
     out->length = in->length;
     out->target = in->target;
diff --git a/libdedx/dedx_mstar.h b/libdedx/dedx_mstar.h
index c4711a2..8396a61 100644
--- a/libdedx/dedx_mstar.h
+++ b/libdedx/dedx_mstar.h
@@ -1,10 +1,29 @@
-#ifndef DEDX_MSTAR_H_INCLUDED
-#define DEDX_MSTAR_H_INCLUDED
+#ifndef DEDX_MSTAR_H
+#define DEDX_MSTAR_H
+
+#include "dedx.h"
 #include "dedx_mpaul.h"
 #include "dedx_stopping_data.h"
 
-void _evaluate_compound_state_mstar(dedx_config *config, int *err);
-void _dedx_convert_energy_to_mstar(
-    stopping_data *in, stopping_data *out, char state, dedx_config *config, float *energy);
+/** @brief Resolve MSTAR auto modes `a`/`b` to the actual condensed or gaseous mode.
+ *  @param[in,out] config  Configuration whose `mstar_mode` may be rewritten.
+ *  @param[out]    err     Error code; 0 on success.
+ *
+ *  MSTAR mode semantics follow the original `MSTAR1` routine:
+ *  `a` selects `c` for condensed targets and `g` for gaseous targets,
+ *  while `b` selects `d` for condensed targets and `h` for gaseous targets.
+ */
+void dedx_internal_evaluate_compound_state_mstar(dedx_config *config, int *err);
+
+/** @brief Convert alpha-particle stopping data to heavy-ion MSTAR stopping data.
+ *  @param[in]     in      Input alpha stopping-power table.
+ *  @param[out]    out     Output heavy-ion stopping-power table.
+ *  @param[in]     state   Requested MSTAR mode (`a`,`b`,`c`,`d`,`g`,`h`).
+ *  @param[in]     config  Resolved configuration, including target state when available.
+ *  @param[in,out] energy  Energy grid, converted in place from alpha to ion scaling.
+ *  @param[out]    err     Error code; 0 on success.
+ */
+void dedx_internal_convert_energy_to_mstar(
+    stopping_data *in, stopping_data *out, char state, dedx_config *config, float *energy, int *err);
 
-#endif // DEDX_MSTAR_H_INCLUDED
+#endif // DEDX_MSTAR_H
diff --git a/libdedx/dedx_periodic_table.c b/libdedx/dedx_periodic_table.c
index ef42600..fa7ca4c 100644
--- a/libdedx/dedx_periodic_table.c
+++ b/libdedx/dedx_periodic_table.c
@@ -1,6 +1,7 @@
 #include "dedx_periodic_table.h"
 
-float _dedx_get_atom_charge(int id, int *err) {
+float dedx_internal_get_atom_charge(int id, int *err) {
+    *err = DEDX_OK;
     if (id < 113) {
         return id;
     }
@@ -8,7 +9,7 @@ float _dedx_get_atom_charge(int id, int *err) {
     return -1;                      /* LCOV_EXCL_LINE */
 }
 
-float _dedx_get_atom_mass(int id, int *err) {
+float dedx_internal_get_atom_mass(int id, int *err) {
     *err = DEDX_OK;
     if (id < 113)
         return dedx_amu[id - 1];
@@ -16,7 +17,7 @@ float _dedx_get_atom_mass(int id, int *err) {
     return -1;                      /* LCOV_EXCL_LINE */
 }
 
-int _dedx_get_nucleon(int id, int *err) {
+int dedx_internal_get_nucleon(int id, int *err) {
     *err = DEDX_OK;
     if (id < 113)
         return dedx_nucl[id - 1];
diff --git a/libdedx/dedx_periodic_table.h b/libdedx/dedx_periodic_table.h
index 2c654a9..f7d3be7 100644
--- a/libdedx/dedx_periodic_table.h
+++ b/libdedx/dedx_periodic_table.h
@@ -1,5 +1,5 @@
-#ifndef DEDX_PERIODIC_TABLE_H_INCLUDED
-#define DEDX_PERIODIC_TABLE_H_INCLUDED
+#ifndef DEDX_PERIODIC_TABLE_H
+#define DEDX_PERIODIC_TABLE_H
 
 #include "dedx_error.h"
 
@@ -17,9 +17,9 @@ static const float dedx_amu[112] = {
     258,         259,       262,        267,      268,         271,     272,       270,     276,        281,
     280,         285};
 
-float _dedx_get_atom_charge(int id, int *err);
+float dedx_internal_get_atom_charge(int id, int *err);
 
-float _dedx_get_atom_mass(int id, int *err);
+float dedx_internal_get_atom_mass(int id, int *err);
 
 static const int dedx_nucl[112] = {
     1,   4,   7,   9,   11,  12,  14,  16,  19,  20,  23,  24,  27,  28,  31,  32,  35,  40,  39,  40,  45,  48,  51,
@@ -28,6 +28,6 @@ static const int dedx_nucl[112] = {
     174, 175, 180, 181, 184, 185, 192, 193, 195, 197, 202, 205, 208, 209, 209, 210, 222, 223, 226, 227, 232, 231, 238,
     237, 244, 243, 247, 247, 251, 252, 257, 258, 259, 262, 267, 268, 271, 272, 270, 276, 281, 280, 285};
 
-int _dedx_get_nucleon(int id, int *err);
+int dedx_internal_get_nucleon(int id, int *err);
 
-#endif // DEDX_PERIODIC_TABLE_H_INCLUDED
+#endif // DEDX_PERIODIC_TABLE_H
diff --git a/libdedx/dedx_program_const.h b/libdedx/dedx_program_const.h
index 7a3abac..6548172 100644
--- a/libdedx/dedx_program_const.h
+++ b/libdedx/dedx_program_const.h
@@ -1,10 +1,9 @@
-#ifndef DEDX_PROGRAM_CONST_H_INCLUDED
-#define DEDX_PROGRAM_CONST_H_INCLUDED
+#ifndef DEDX_PROGRAM_CONST_H
+#define DEDX_PROGRAM_CONST_H
 
-/* ok, this is a bit overkill... */
+/* Internal lookup tables for program names, versions, and supported IDs. */
 
-/* Hey Jakob, whats the purpose with  "" placeholder ? */
-const char dedx_program_table[110][20] = {
+static const char dedx_program_table[110][20] = {
     // clang-format off
   "(N/A)","ASTAR","PSTAR","ESTAR","MSTAR","ICRU73_OLD","ICRU73","ICRU49","","ICRU", // 0 - 9
   "","","","","","","","","","",  // 10-19
@@ -18,7 +17,7 @@ const char dedx_program_table[110][20] = {
   "","","","","","","","","","",  // 90-99
   "DEFAULT","BETHE_EXT00","","","","","","","",""   // 100-109  // clang-format on
 };
-const char dedx_program_version_table[110][20] = {  // clang-format off
+static const char dedx_program_version_table[110][20] = {  // clang-format off
   "(N/A)","1.2","1.2","ESTAR","3.12","2005/06","2005/06-REV","ICRU Report 49 D","","1.0", // 0 - 9
   "","","","","","","","","","",  // 10-19
   "","","","","","","","","","",  // 20-29
@@ -32,14 +31,14 @@ const char dedx_program_version_table[110][20] = {  // clang-format off
   "1.0","","","","","","","","",""   // 100-109  // clang-format on
 };
 
-const int dedx_available_programs[20] = {
+static const int dedx_available_programs[20] = {
   DEDX_ASTAR, DEDX_PSTAR,
   DEDX_ESTAR, DEDX_MSTAR,
   DEDX_ICRU73_OLD, DEDX_ICRU73,
   DEDX_ICRU49, DEDX_ICRU, DEDX_DEFAULT,
   DEDX_BETHE_EXT00, -1};
 
-const int dedx_program_available_ions[110][20] = {
+static const int dedx_program_available_ions[110][20] = {
   // currently excluding BETHE, which is handled in dedx.c
   {1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,-1}, // all ions
   {2,-1}, // ASTAR
@@ -49,10 +48,10 @@ const int dedx_program_available_ions[110][20] = {
   {3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,-1}, // ICRU73_OLD
   {3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,-1}, // ICRU73
   {1,2,-1}, // ICRU49
-  {-1}, // Jakobs magic placeholder
+  {-1}, // reserved internal slot
   {1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,-1} // ICRU
 };
-const int dedx_program_available_materials[110][290] = {
+static const int dedx_program_available_materials[110][290] = {
   // currently excluding BETHE, which is handled in dedx.c
   {1,2,3,4,5,6,7,8,9,
    10,11,12,13,14,15,16,17,18,19,
@@ -121,14 +120,14 @@ const int dedx_program_available_materials[110][290] = {
    73,74,78,79,82,92,99,101,103,104,106,111,119,120,126,130,134,138,139,
    141,155,160,169,179,185,189,191,197,200,201,202,203,204,209,213,215,216,219,
    221,222,223,225,226,227,232,238,245,252,255,263,264,266,276,277,906,-1},
-  {-1}, // Jakobs magic placeholder
+  {-1}, // reserved internal slot
   {1,2,4,5,6,7,8,10,13,14,18,22,26,28,29,32,36,40,42,47,50,54,64, // ICRU
    73,74,78,79,82,92,99,101,103,104,106,111,119,120,126,130,134,138,139,
    141,155,160,169,179,185,189,191,197,200,201,202,203,204,209,213,215,216,219,
    221,222,223,225,226,227,232,238,245,252,255,263,264,266,276,277,906,-1}
 };
 
-const char dedx_material_table[1000][40] = {
+static const char dedx_material_table[1000][40] = {
   "(N/A)",
   "HYDROGEN", "HELIUM", "LITHIUM", "BERYLLIUM", "BORON", "CARBON",
   "NITROGEN", "OXYGEN", "FLUORINE", "NEON", "SODIUM", "MAGNESIUM",
@@ -242,7 +241,7 @@ const char dedx_material_table[1000][40] = {
   "","","","","","","","","",""  // last element is 999
 };
 
-const char dedx_ion_table[120][40] = {
+static const char dedx_ion_table[120][40] = {
   "(N/A)",
   "HYDROGEN", "HELIUM", "LITHIUM", "BERYLLIUM", "BORON", "CARBON",
   "NITROGEN", "OXYGEN", "FLUORINE", "NEON", "SODIUM", "MAGNESIUM",
@@ -265,4 +264,4 @@ const char dedx_ion_table[120][40] = {
   "ROENTGENIUM", "COPERNICUM" // Copernicum = 112
 };
 
-#endif // DEDX_PROGRAM_CONST_H_INCLUDED
+#endif // DEDX_PROGRAM_CONST_H
diff --git a/libdedx/dedx_spline.c b/libdedx/dedx_spline.c
index 98bc31f..8e2a9d6 100644
--- a/libdedx/dedx_spline.c
+++ b/libdedx/dedx_spline.c
@@ -16,17 +16,10 @@
 */
 #include "dedx_spline.h"
 
-int _dedx_linear_search(_dedx_spline_base *coef, float value) {
-    int i;
-    for (i = 0; i < _DEDX_MAXELEMENTS; i++) {
-        if (coef[i].x >= value) {
-            break;
-        }
-    }
-    return i;
-}
+#include <math.h>
+#include <stddef.h>
 
-int _dedx_binary_search(_dedx_spline_base *coef, float value, int n) {
+static int binary_search(_dedx_spline_base *coef, float value, int n) {
     int head = n - 1;
     int tail = 0;
     int guess = n / 2;
@@ -41,13 +34,16 @@ int _dedx_binary_search(_dedx_spline_base *coef, float value, int n) {
     return guess;
 }
 
-void _dedx_calculate_coefficients(_dedx_spline_base *coef, float *energy, float *stopping, int n) {
+void dedx_internal_calculate_coefficients(_dedx_spline_base *coef, float *energy, float *stopping, int n) {
     int i;
-    float h[_DEDX_MAXELEMENTS];
-    float alpha[_DEDX_MAXELEMENTS];
-    float l[_DEDX_MAXELEMENTS];
-    float my[_DEDX_MAXELEMENTS];
-    float z[_DEDX_MAXELEMENTS];
+    float h[DEDX_MAX_ELEMENTS];
+    float alpha[DEDX_MAX_ELEMENTS];
+    float l[DEDX_MAX_ELEMENTS];
+    float my[DEDX_MAX_ELEMENTS];
+    float z[DEDX_MAX_ELEMENTS];
+
+    if (n < 2)
+        return;
 
     l[0] = 1;
     my[0] = 0;
@@ -79,20 +75,17 @@ void _dedx_calculate_coefficients(_dedx_spline_base *coef, float *energy, float
     }
 }
 
-float _dedx_evaluate_spline(_dedx_spline_base *coef, float x, _dedx_lookup_accelerator *acc, int n) {
+float dedx_internal_evaluate_spline(_dedx_spline_base *coef, float x, _dedx_lookup_accelerator *acc, int n) {
     int i;
     int lookup = 1;
     if (acc != NULL) {
-        // Next line looks weird.
-        // if((coef[acc->cache].x <= x) & x < coef[acc->cache+1].x)
-        // TBC: Just a guess what Jakob meant?
         if ((coef[acc->cache].x <= x) && (x < coef[acc->cache + 1].x)) {
             lookup = 0;
             i = acc->cache;
         }
     }
     if (lookup) {
-        i = _dedx_binary_search(coef, x, n);
+        i = binary_search(coef, x, n);
         if (acc != NULL)
             acc->cache = i;
     }
diff --git a/libdedx/dedx_spline.h b/libdedx/dedx_spline.h
index 73acc47..d346040 100644
--- a/libdedx/dedx_spline.h
+++ b/libdedx/dedx_spline.h
@@ -1,14 +1,24 @@
-#ifndef DEDX_SPLINE_H_INCLUDED
-#define DEDX_SPLINE_H_INCLUDED
-#include <math.h>
-#include <stdio.h>
-#include <stdlib.h>
+#ifndef DEDX_SPLINE_H
+#define DEDX_SPLINE_H
 
 #include "dedx_lookup_accelerator.h"
 #include "dedx_spline_base.h"
 
-void _dedx_calculate_coefficients(_dedx_spline_base *coef, float *energy, float *stopping, int n);
+/** @brief Build cubic-spline coefficients for one stopping-power table.
+ *  @param[out] coef      Spline coefficient array to populate.
+ *  @param[in]  energy    Monotonic energy grid.
+ *  @param[in]  stopping  Stopping-power values on @p energy.
+ *  @param[in]  n         Number of valid grid points.
+ */
+void dedx_internal_calculate_coefficients(_dedx_spline_base *coef, float *energy, float *stopping, int n);
 
-float _dedx_evaluate_spline(_dedx_spline_base *coef, float x, _dedx_lookup_accelerator *acc, int n);
+/** @brief Evaluate a precomputed cubic spline at one energy value.
+ *  @param[in]      coef  Spline coefficients.
+ *  @param[in]      x     Energy at which to evaluate the spline.
+ *  @param[in,out]  acc   Optional lookup cache for repeated evaluations; may be NULL.
+ *  @param[in]      n     Number of valid spline intervals.
+ *  @return Interpolated stopping power at @p x.
+ */
+float dedx_internal_evaluate_spline(_dedx_spline_base *coef, float x, _dedx_lookup_accelerator *acc, int n);
 
-#endif // DEDX_SPLINE_H_INCLUDED
+#endif // DEDX_SPLINE_H
diff --git a/libdedx/dedx_spline_base.h b/libdedx/dedx_spline_base.h
index b7dccf6..137595e 100644
--- a/libdedx/dedx_spline_base.h
+++ b/libdedx/dedx_spline_base.h
@@ -1,6 +1,14 @@
-#ifndef DEDX_SPLINE_BASE_H_INCLUDED
-#define DEDX_SPLINE_BASE_H_INCLUDED
+#ifndef DEDX_SPLINE_BASE_H
+#define DEDX_SPLINE_BASE_H
 
-#include "dedx.h"
+#include "dedx_elements.h"
 
-#endif // DEDX_SPLINE_BASE_H_INCLUDED
+typedef struct {
+    float a;
+    float b;
+    float c;
+    float d;
+    float x;
+} _dedx_spline_base;
+
+#endif // DEDX_SPLINE_BASE_H
diff --git a/libdedx/dedx_split.c b/libdedx/dedx_split.c
index 658bd88..09c0eae 100644
--- a/libdedx/dedx_split.c
+++ b/libdedx/dedx_split.c
@@ -18,21 +18,27 @@
 
 #include <stdlib.h>
 
-char **_dedx_split(char *string, char splitCaracter, unsigned int *items, unsigned int strLength) {
+char **dedx_internal_split(char *string, char split_character, unsigned int *items, unsigned int str_length) {
     int i = 0;
     int j = 0;
-    char **words;
-    words = malloc(10 * sizeof(char *));
+    char **words = (char **) malloc(10 * sizeof(char *));
     *items = 0;
     char word[20] = "";
 
+    if (words == NULL)
+        return NULL;
+
     for (i = 0; i < 10; i++) {
-        words[i] = (char *) malloc(20 * sizeof(char));
+        words[i] = malloc(20 * sizeof(char));
+        if (words[i] == NULL) {
+            dedx_internal_free_split_temp(words, i);
+            return NULL;
+        }
     }
 
-    for (i = 0; i < strLength; i++) {
+    for (i = 0; i < str_length; i++) {
 
-        if (string[i] == splitCaracter) {
+        if (string[i] == split_character) {
 
             for (j = 0; j < 20; j++) {
                 words[*items][j] = word[j];
@@ -55,9 +61,11 @@ char **_dedx_split(char *string, char splitCaracter, unsigned int *items, unsign
     return words;
 }
 
-void _dedx_free_split_temp(char **temp) {
-    int i = 0;
-    for (i = 0; i < 10; i++)
+void dedx_internal_free_split_temp(char **temp, unsigned int count) {
+    unsigned int i = 0;
+    if (temp == NULL)
+        return;
+    for (i = 0; i < count; i++)
         free(temp[i]);
-    free(temp);
+    free((void *) temp);
 }
diff --git a/libdedx/dedx_split.h b/libdedx/dedx_split.h
index d2aaef4..1a813ab 100644
--- a/libdedx/dedx_split.h
+++ b/libdedx/dedx_split.h
@@ -1,6 +1,22 @@
-#ifndef DEDX_SPLIT_H_INCLUDED
-#define DEDX_SPLIT_H_INCLUDED
-char **_dedx_split(char *string, char splitCaracter, unsigned int *items, unsigned int strLength);
-void _dedx_free_split_temp(char **temp);
+#ifndef DEDX_SPLIT_H
+#define DEDX_SPLIT_H
 
-#endif // DEDX_SPLIT_H_INCLUDED
+/** @brief Split a fixed-width line buffer into heap-allocated tokens.
+ *  @param[in]  string           Input buffer to split.
+ *  @param[in]  split_character  Delimiter character.
+ *  @param[out] items            Number of parsed tokens.
+ *  @param[in]  str_length       Maximum number of characters to scan in @p string.
+ *  @return Array of 10 heap-allocated token buffers, or NULL on allocation failure.
+ *
+ *  Each token buffer has room for up to 20 characters including the terminator.
+ *  The returned array must be released with dedx_internal_free_split_temp().
+ */
+char **dedx_internal_split(char *string, char split_character, unsigned int *items, unsigned int str_length);
+
+/** @brief Free a token array returned by dedx_internal_split().
+ *  @param[in] temp   Token array to free; NULL is allowed.
+ *  @param[in] count  Number of allocated token buffers in @p temp.
+ */
+void dedx_internal_free_split_temp(char **temp, unsigned int count);
+
+#endif // DEDX_SPLIT_H
diff --git a/libdedx/dedx_stopping_data.h b/libdedx/dedx_stopping_data.h
index eee5e1f..c97fa6f 100644
--- a/libdedx/dedx_stopping_data.h
+++ b/libdedx/dedx_stopping_data.h
@@ -1,13 +1,13 @@
-#ifndef DEDX_STOPPING_DATA_H_INCLUDED
-#define DEDX_STOPPING_DATA_H_INCLUDED
+#ifndef DEDX_STOPPING_DATA_H
+#define DEDX_STOPPING_DATA_H
 
-#include "dedx.h"
+#include "dedx_elements.h"
 
 typedef struct {
     int target;
     int ion;
     unsigned int length;
-    float data[_DEDX_MAXELEMENTS];
+    float data[DEDX_MAX_ELEMENTS];
 } stopping_data;
 
-#endif // DEDX_STOPPING_DATA_H_INCLUDED
+#endif // DEDX_STOPPING_DATA_H
diff --git a/libdedx/dedx_tools.c b/libdedx/dedx_tools.c
index 3b3d785..ec153fe 100644
--- a/libdedx/dedx_tools.c
+++ b/libdedx/dedx_tools.c
@@ -7,21 +7,19 @@
 #include "dedx_file_access.h"
 
 typedef struct {
-    // int id;
-    //   double A;
     dedx_workspace *ws;
     dedx_config *cfg;
-} _dedx_tools_settings;
-
-double _dedx_adapt24(double (*func)(double x, _dedx_tools_settings *set),
-                     _dedx_tools_settings *set,
-                     double a,
-                     double b,
-                     double f2,
-                     double f3,
-                     double acc,
-                     double eps,
-                     double *err) {
+} dedx_tools_settings;
+
+static double adapt24(double (*func)(double x, dedx_tools_settings *set),
+                      dedx_tools_settings *set,
+                      double a,
+                      double b,
+                      double f2,
+                      double f3,
+                      double acc,
+                      double eps,
+                      double *err) {
     double h = b - a;
     double f1 = (*func)(a + h / 6, set);
     double f4 = (*func)(a + h * 5 / 6, set);
@@ -31,44 +29,45 @@ double _dedx_adapt24(double (*func)(double x, _dedx_tools_settings *set),
     *err = fabs(q4 - q2) / 3;
     if (*err < tol)
         return q4;
-    else {
-        acc /= sqrt(2);
-        double mid = (a + b) / 2;
-        double el = 0;
-        double er = 0;
-        double ql = _dedx_adapt24(func, set, a, mid, f1, f2, acc, eps, &el);
-        double qr = _dedx_adapt24(func, set, mid, b, f3, f4, acc, eps, &er);
-        *err = sqrt(el * el + er * er);
-        return ql + qr;
-    }
+
+    acc /= sqrt(2);
+    double mid = (a + b) / 2;
+    double el = 0;
+    double er = 0;
+    double ql = adapt24(func, set, a, mid, f1, f2, acc, eps, &el);
+    double qr = adapt24(func, set, mid, b, f3, f4, acc, eps, &er);
+    *err = sqrt(el * el + er * er);
+    return ql + qr;
 }
 
-double _dedx_adapt(double (*func)(double x, _dedx_tools_settings *set),
-                   _dedx_tools_settings *set,
-                   double a,
-                   double b,
-                   double acc,
-                   double eps,
-                   double *err) {
+static double adapt(double (*func)(double x, dedx_tools_settings *set),
+                    dedx_tools_settings *set,
+                    double a,
+                    double b,
+                    double acc,
+                    double eps,
+                    double *err) {
     double h = b - a;
-    return _dedx_adapt24(func, set, a, b, (*func)(a + h / 3, set), (*func)(a + 2 * h / 3, set), acc, eps, err);
+    return adapt24(func, set, a, b, (*func)(a + h / 3, set), (*func)(a + 2 * h / 3, set), acc, eps, err);
 }
 
-double _dedx_adapt_stp(double energy, _dedx_tools_settings *set) {
+static double adapt_stp(double energy, dedx_tools_settings *set) {
     int err = 0;
-    // TODO: i am sure sure if its cfg->ion_a or ion_A !!
-    return 1 / dedx_get_stp(set->ws, set->cfg, energy / (set->cfg->ion_a), &err);
+    double stp = dedx_get_stp(set->ws, set->cfg, energy / (set->cfg->ion_a), &err);
+    if (err != 0 || stp == 0.0)
+        return INFINITY;
+    return 1.0 / stp;
 }
 
-double _dedx_find_min_stp_func(double x, _dedx_tools_settings *set) {
+static double find_min_stp_func(double x, dedx_tools_settings *set) {
     int err = 0;
-    float stp = 1 / dedx_get_stp(set->ws, set->cfg, x, &err);
-    if (err != 0)
+    double stp = dedx_get_stp(set->ws, set->cfg, x, &err);
+    if (err != 0 || stp == 0.0)
         return INFINITY;
-    return stp;
+    return 1.0 / stp;
 }
 
-double _dedx_find_min(double (*func)(double x, _dedx_tools_settings *set), _dedx_tools_settings *set, double acc) {
+static double find_min(double (*func)(double x, dedx_tools_settings *set), dedx_tools_settings *set, double acc) {
     double x[] = {0.01, 10};
     double f[] = {func(x[0], set), func(x[1], set)};
     int i = 0;
@@ -131,17 +130,16 @@ double dedx_get_inverse_stp(dedx_workspace *ws, dedx_config *config, float stp,
         return -1;
     }
     double acc = 1e-5;
-    _dedx_tools_settings set;
+    dedx_tools_settings set;
     set.ws = ws;
     if (*err != 0)
         return -1;
     dedx_load_config(ws, config, err);
-    //	set.id = config->cfg_id;
     set.cfg = config;
 
     if (*err != 0)
         return -1;
-    double max = _dedx_find_min(_dedx_find_min_stp_func, &set, acc * 100);
+    double max = find_min(find_min_stp_func, &set, acc * 100);
     double x1;
     double x2;
     double x_temp;
@@ -163,7 +161,6 @@ double dedx_get_inverse_stp(dedx_workspace *ws, dedx_config *config, float stp,
             x1 = x_temp;
         }
     }
-    // free(cfg); might be used later again
     return (x1 + x2) / 2;
 }
 
@@ -175,9 +172,9 @@ double dedx_get_csda(dedx_workspace *ws, dedx_config *config, float energy, int
     double calculation_error = 0;
     double acc = 1e-6;
     double eps = 1e-6;
-    _dedx_tools_settings set;
+    dedx_tools_settings set;
     double range = 0.0;
-    double A = config->ion_a; // TODO: or nucleon number?
+    double A = config->ion_a;
 
     if (*err != 0)
         return -1;
@@ -186,20 +183,19 @@ double dedx_get_csda(dedx_workspace *ws, dedx_config *config, float energy, int
     if (*err != 0)
         return -1;
     set.cfg = config;
-    // set.A = config->nucleon_number;
     set.ws = ws;
-    range = _dedx_adapt(_dedx_adapt_stp,
-                        &set,
-                        dedx_get_min_energy(config->program, config->ion) * A,
-                        energy * A,
-                        acc,
-                        eps,
-                        &calculation_error);
+    range = adapt(adapt_stp,
+                  &set,
+                  dedx_get_min_energy(config->program, config->ion) * A,
+                  energy * A,
+                  acc,
+                  eps,
+                  &calculation_error);
     return range;
 }
 
-float _conversion_factor(const int old_unit, const int new_unit, const int material, int *err) {
-    const float density = _dedx_read_density(material, err);
+static float conversion_factor(const int old_unit, const int new_unit, const int material, int *err) {
+    const float density = dedx_internal_read_density(material, err);
 
     float conversion_rate;
 
@@ -242,11 +238,12 @@ int convert_units(const int old_unit,
                   const float *old_values,
                   float *new_values) {
     int err = 0;
+    int i;
     if (old_unit == new_unit)
         return err;
 
-    float conversion_rate = _conversion_factor(old_unit, new_unit, material, &err);
-    for (int i = 0; i < no_of_points; i++) {
+    float conversion_rate = conversion_factor(old_unit, new_unit, material, &err);
+    for (i = 0; i < no_of_points; i++) {
         new_values[i] = old_values[i] * conversion_rate;
     }
     return err;
diff --git a/libdedx/dedx_tools.h b/libdedx/dedx_tools.h
index bdd7ae4..15df564 100644
--- a/libdedx/dedx_tools.h
+++ b/libdedx/dedx_tools.h
@@ -1,5 +1,5 @@
-#ifndef DEDX_TOOLS_H_INCLUDED
-#define DEDX_TOOLS_H_INCLUDED
+#ifndef DEDX_TOOLS_H
+#define DEDX_TOOLS_H
 
 /**
  * @file dedx_tools.h
@@ -54,10 +54,6 @@ double dedx_get_inverse_stp(dedx_workspace *ws, dedx_config *config, float stp,
  */
 double dedx_get_inverse_csda(dedx_workspace *ws, dedx_config *config, float range, int *err);
 
-/** @cond INTERNAL */
-float _conversion_factor(const int old_unit, const int new_unit, const int material, int *err);
-/** @endcond */
-
 /** @brief Convert an array of stopping power values between unit systems.
  *
  *  @param[in]  old_unit      Source unit (dedx_stp_units).
@@ -75,4 +71,4 @@ int convert_units(const int old_unit,
                   const float *old_values,
                   float *new_values);
 
-#endif // DEDX_TOOLS_H_INCLUDED
+#endif // DEDX_TOOLS_H
diff --git a/libdedx/dedx_validate.c b/libdedx/dedx_validate.c
index ab7a7a1..bd81d86 100644
--- a/libdedx/dedx_validate.c
+++ b/libdedx/dedx_validate.c
@@ -1,77 +1,113 @@
 #include "dedx_validate.h"
 
 #include <math.h>
+#include <stdlib.h>
 
 #include "dedx_file_access.h"
 #include "dedx_periodic_table.h"
 
-int _dedx_set_names(dedx_config *config, int *err) {
+int dedx_internal_set_names(dedx_config *config, int *err) {
     if (config->target != 0)
         config->target_name = dedx_get_material_name(config->target);
 
     config->ion_name = dedx_get_ion_name(config->ion);
     config->program_name = dedx_get_program_name(config->program);
+    *err = DEDX_OK;
     return 0;
 }
 
-int _dedx_validate_rho(dedx_config *config, int *err) {
+int dedx_internal_validate_rho(dedx_config *config, int *err) {
     if (config->rho <= 0.0 && config->target != 0) {
-        config->rho = _dedx_read_density(config->target, err);
+        config->rho = dedx_internal_read_density(config->target, err);
     } else if (config->rho <= 0.0 && config->target == 0 && config->program >= 100) {
         *err = DEDX_ERR_RHO_REQUIRED;
     }
     return 0;
 }
 
-int _dedx_evaluate_i_pot(dedx_config *config, int *err) {
+int dedx_internal_evaluate_i_pot(dedx_config *config, int *err) {
     int i;
+
     if (config->elements_i_value == NULL && config->target != 0) {
         if (config->i_value == 0.0) {
-            config->i_value = _dedx_get_i_value(config->target, config->compound_state, err);
+            config->i_value = dedx_internal_get_i_value(config->target, config->compound_state, err);
         }
         if (*err != 0)
             return -1;
     } else if (config->i_value == 0.0 && config->target == 0 && config->elements_i_value == NULL) {
+        if (config->elements_length == 0) {
+            *err = DEDX_ERR_TARGET_NOT_FOUND;
+            return -1;
+        }
         config->elements_i_value = calloc(config->elements_length, sizeof(float));
+        if (config->elements_i_value == NULL) {
+            *err = DEDX_ERR_NO_MEMORY;
+            return -1;
+        }
         for (i = 0; i < config->elements_length; i++) {
-            config->elements_i_value[i] = _dedx_get_i_value(config->elements_id[i], config->compound_state, err);
+            config->elements_i_value[i] =
+                dedx_internal_get_i_value(config->elements_id[i], config->compound_state, err);
+            if (*err != 0)
+                return -1;
         }
     }
     if (config->elements_id != NULL && config->elements_i_value == NULL) {
-        _dedx_calculate_element_i_pot(config, err);
+        dedx_internal_calculate_element_i_pot(config, err);
+        if (*err != 0)
+            return -1;
     }
     if (config->elements_i_value != NULL && config->i_value == 0.0) {
+        if (config->elements_id == NULL || config->elements_mass_fraction == NULL) {
+            *err = DEDX_ERR_INCONSISTENT_COMPOUND;
+            return -1;
+        }
         float charge_avg = 0.0;
         for (i = 0; i < config->elements_length; i++) {
             config->i_value += config->elements_mass_fraction[i] * log(config->elements_i_value[i])
-                               * config->elements_id[i] / _dedx_get_atom_mass(config->elements_id[i], err);
+                               * config->elements_id[i] / dedx_internal_get_atom_mass(config->elements_id[i], err);
             charge_avg += config->elements_mass_fraction[i] * config->elements_id[i]
-                          / _dedx_get_atom_mass(config->elements_id[i], err);
+                          / dedx_internal_get_atom_mass(config->elements_id[i], err);
+            if (*err != 0)
+                return -1;
         }
         config->i_value = exp(config->i_value / charge_avg);
     }
     return 0;
 }
 
-int _dedx_evaluate_compound(dedx_config *config, int *err) {
+int dedx_internal_evaluate_compound(dedx_config *config, int *err) {
     int i = 0;
+
     if (config->target > 0 && config->target <= 99) {
+        *err = DEDX_OK;
         return 0;
     }
+
     config->bragg_used = 1;
-    //_evaluate_compound_state_mstar(config,err);
     if (config->elements_id == NULL) {
         unsigned int compos_len;
         float composition[20][2];
-        _dedx_get_composition(config->target, composition, &compos_len, err);
-        if (*err != 0)
+
+        dedx_internal_get_composition(config->target, composition, &compos_len, err);
+        if (*err != 0) {
             return -1;
+        }
         if (compos_len == 0) { /* LCOV_EXCL_START */
             *err = DEDX_ERR_TARGET_NOT_FOUND;
             return -1;
         } /* LCOV_EXCL_STOP */
+
         config->elements_id = (int *) malloc(sizeof(int) * compos_len);
         config->elements_mass_fraction = (float *) malloc(sizeof(float) * compos_len);
+        if (config->elements_id == NULL || config->elements_mass_fraction == NULL) {
+            free(config->elements_id);
+            free(config->elements_mass_fraction);
+            config->elements_id = NULL;
+            config->elements_mass_fraction = NULL;
+            *err = DEDX_ERR_NO_MEMORY;
+            return -1;
+        }
+
         for (i = 0; i < compos_len; i++) {
             config->elements_id[i] = (int) composition[i][0];
             config->elements_mass_fraction[i] = composition[i][1];
@@ -79,57 +115,76 @@ int _dedx_evaluate_compound(dedx_config *config, int *err) {
         config->elements_length = compos_len;
     } else if (config->elements_mass_fraction == NULL && config->elements_atoms != NULL) {
         int length = config->elements_length;
-        int *compos = config->elements_atoms;
+        int *atoms_per_element = config->elements_atoms;
         float *density = malloc(sizeof(float) * length);
         float *weight = malloc(sizeof(float) * length);
         float f, sum = 0;
+
+        if (length == 0) {
+            free(density);
+            free(weight);
+            *err = DEDX_ERR_TARGET_NOT_FOUND;
+            return -1;
+        }
+        if (density == NULL || weight == NULL) {
+            free(density);
+            free(weight);
+            *err = DEDX_ERR_NO_MEMORY;
+            return -1;
+        }
         for (i = 0; i < length; i++) {
-            f = _dedx_get_atom_mass(config->elements_id[i], err);
+            f = dedx_internal_get_atom_mass(config->elements_id[i], err);
             if (*err != 0) {
                 free(density);
                 free(weight);
                 return -1;
             }
             density[i] = f;
-            sum += compos[i] * f;
+            sum += atoms_per_element[i] * f;
         }
 
         for (i = 0; i < length; i++) {
-            weight[i] = compos[i] * density[i] / sum;
+            weight[i] = atoms_per_element[i] * density[i] / sum;
         }
         free(density);
         config->elements_mass_fraction = weight;
     } else {
+        *err = DEDX_ERR_INCONSISTENT_COMPOUND;
         return -1;
     }
     return 0;
 }
 
-int _dedx_validate_config(dedx_config *config, int *err) {
-    _dedx_validate_rho(config, err);
-    if (*err != 0)
+int dedx_internal_validate_config(dedx_config *config, int *err) {
+    dedx_internal_validate_rho(config, err);
+    if (*err != 0) {
         return -1;
+    }
+
     if (config->program >= 100) {
         // Order is important
-        _dedx_evaluate_compound(config, err);
+        dedx_internal_evaluate_compound(config, err);
         if (*err != 0)
             return -1;
-        _dedx_validate_state(config, err);
+        dedx_internal_validate_state(config, err);
         if (*err != 0)
             return -1;
-        _dedx_evaluate_i_pot(config, err);
+        dedx_internal_evaluate_i_pot(config, err);
         if (*err != 0)
             return -1;
     }
+
     if (config->target == 0 && config->elements_id != NULL) {
-        _dedx_evaluate_compound(config, err);
+        dedx_internal_evaluate_compound(config, err);
+        if (*err != 0)
+            return -1;
     }
     return 0;
 }
 
-int _dedx_validate_state(dedx_config *config, int *err) {
+int dedx_internal_validate_state(dedx_config *config, int *err) {
     if (config->compound_state == DEDX_DEFAULT_STATE) {
-        if (_dedx_target_is_gas(config->target, err)) {
+        if (dedx_internal_target_is_gas(config->target, err)) {
             config->compound_state = DEDX_GAS;
         } else {
             config->compound_state = DEDX_CONDENSED;
@@ -138,24 +193,42 @@ int _dedx_validate_state(dedx_config *config, int *err) {
     return 0;
 }
 
-int _dedx_calculate_element_i_pot(dedx_config *config, int *err) {
+int dedx_internal_calculate_element_i_pot(dedx_config *config, int *err) {
     int i;
     float charge_avg = 0;
     float avg_pot = 0;
     float log_x, i_pot_x;
     int target;
+
+    if (config->elements_length == 0) {
+        *err = DEDX_ERR_TARGET_NOT_FOUND;
+        return -1;
+    }
+
     for (i = 0; i < config->elements_length; i++) {
         target = config->elements_id[i];
-        charge_avg += config->elements_mass_fraction[i] * target / _dedx_get_atom_mass(target, err);
-        avg_pot += config->elements_mass_fraction[i] * target / _dedx_get_atom_mass(target, err)
-                   * log(_dedx_get_i_value(target, config->compound_state, err));
+        charge_avg += config->elements_mass_fraction[i] * target / dedx_internal_get_atom_mass(target, err);
+        avg_pot += config->elements_mass_fraction[i] * target / dedx_internal_get_atom_mass(target, err)
+                   * log(dedx_internal_get_i_value(target, config->compound_state, err));
+        if (*err != 0)
+            return -1;
     }
+
     log_x = log(config->i_value);
     log_x -= avg_pot / charge_avg;
     i_pot_x = exp(log_x);
     config->elements_i_value = (float *) malloc(sizeof(float) * config->elements_length);
+    if (config->elements_i_value == NULL) {
+        *err = DEDX_ERR_NO_MEMORY;
+        return -1;
+    }
+
     for (i = 0; i < config->elements_length; i++) {
-        config->elements_i_value[i] = _dedx_get_i_value(config->elements_id[i], config->compound_state, err) * i_pot_x;
+        config->elements_i_value[i] =
+            dedx_internal_get_i_value(config->elements_id[i], config->compound_state, err) * i_pot_x;
+        if (*err != 0)
+            return -1;
     }
+
     return 0;
 }
diff --git a/libdedx/dedx_validate.h b/libdedx/dedx_validate.h
index 41dd211..ce908c4 100644
--- a/libdedx/dedx_validate.h
+++ b/libdedx/dedx_validate.h
@@ -1,13 +1,13 @@
-#ifndef DEDX_VALIDATE_H_INCLUDED
-#define DEDX_VALIDATE_H_INCLUDED
+#ifndef DEDX_VALIDATE_H
+#define DEDX_VALIDATE_H
 #include "dedx.h"
 
-int _dedx_set_names(dedx_config *config, int *err);
-int _dedx_validate_state(dedx_config *config, int *err);
-int _dedx_validate_config(dedx_config *config, int *err);
-int _dedx_evaluate_i_pot(dedx_config *config, int *err);
-int _dedx_validate_rho(dedx_config *config, int *err);
-int _dedx_evaluate_compound(dedx_config *config, int *err);
-int _dedx_calculate_element_i_pot(dedx_config *config, int *err);
+int dedx_internal_set_names(dedx_config *config, int *err);
+int dedx_internal_validate_state(dedx_config *config, int *err);
+int dedx_internal_validate_config(dedx_config *config, int *err);
+int dedx_internal_evaluate_i_pot(dedx_config *config, int *err);
+int dedx_internal_validate_rho(dedx_config *config, int *err);
+int dedx_internal_evaluate_compound(dedx_config *config, int *err);
+int dedx_internal_calculate_element_i_pot(dedx_config *config, int *err);
 
-#endif
+#endif // DEDX_VALIDATE_H
diff --git a/libdedx/dedx_workspace.h b/libdedx/dedx_workspace.h
deleted file mode 100644
index c4d45fe..0000000
--- a/libdedx/dedx_workspace.h
+++ /dev/null
@@ -1,6 +0,0 @@
-#ifndef DEDX_WORKSPACE_H_INCLUDED
-#define DEDX_WORKSPACE_H_INCLUDED
-
-#include "dedx_lookup_data.h"
-
-#endif // DEDX_WORKSPACE_H_INCLUDED
diff --git a/libdedx/dedx_wrappers.c b/libdedx/dedx_wrappers.c
index 01014eb..c993a54 100644
--- a/libdedx/dedx_wrappers.c
+++ b/libdedx/dedx_wrappers.c
@@ -17,9 +17,15 @@
 
 #include "dedx_wrappers.h"
 
+#include <stdlib.h>
+
+#include "dedx_lookup_data.h"
 #include "dedx_periodic_table.h"
 #include "dedx_tools.h"
 
+static dedx_config *allocate_wrapper_config(int program, int ion, int target, int *err);
+static dedx_workspace *allocate_wrapper_workspace(int *err);
+
 void dedx_fill_program_list(int *program_list) {
     /* fill list of available programs, terminated with -1 */
     int i = 0;
@@ -53,23 +59,52 @@ void dedx_fill_ion_list(int program, int *ion_list) {
     ion_list[i] = -1;
 }
 
-// energy - in MeV/nucleon, unless  program is ESTAR, then MeV
-// stp - in MeVcm2/g
+static dedx_config *allocate_wrapper_config(int program, int ion, int target, int *err) {
+    dedx_config *config;
+
+    *err = DEDX_OK;
+    config = (dedx_config *) calloc(1, sizeof(dedx_config));
+    if (config == NULL) {
+        *err = DEDX_ERR_NO_MEMORY;
+        return NULL;
+    }
+
+    config->program = program;
+    config->ion = ion;
+    config->target = target;
+    return config;
+}
+
+static dedx_workspace *allocate_wrapper_workspace(int *err) {
+    return dedx_allocate_workspace(1, err);
+}
+
 int dedx_get_stp_table(
     const int program, const int ion, const int target, const int no_of_points, const float *energies, float *stps) {
     int err = 0;
-    dedx_config *config = (dedx_config *) calloc(1, sizeof(dedx_config));
-    config->target = target;
-    config->ion = ion;
-    config->program = program;
-    dedx_workspace *ws = dedx_allocate_workspace(1, &err);
+    int cleanup_err = DEDX_OK;
+    int i;
+    dedx_config *config = allocate_wrapper_config(program, ion, target, &err);
+    dedx_workspace *ws;
 
     if (err != 0)
         return err;
+    ws = allocate_wrapper_workspace(&err);
+    if (err != 0) {
+        dedx_free_config(config, &cleanup_err);
+        return err;
+    }
     dedx_load_config(ws, config, &err);
+    if (err != 0) {
+        dedx_free_config(config, &cleanup_err);
+        dedx_free_workspace(ws, &cleanup_err);
+        return err;
+    }
 
-    for (int i = 0; i < no_of_points; i++) {
+    for (i = 0; i < no_of_points; i++) {
         stps[i] = dedx_get_stp(ws, config, energies[i], &err);
+        if (err != 0)
+            break;
     }
     dedx_free_config(config, &err);
     dedx_free_workspace(ws, &err);
@@ -78,85 +113,104 @@ int dedx_get_stp_table(
 }
 
 float dedx_get_simple_stp_for_program(const int program, const int ion, const int target, float energy, int *err) {
+    int cleanup_err = DEDX_OK;
     float stp;
-    dedx_config *config = (dedx_config *) calloc(1, sizeof(dedx_config));
-    config->target = target;
-    config->ion = ion;
-    config->program = program;
-    dedx_workspace *ws = dedx_allocate_workspace(1, err);
+    dedx_config *config = allocate_wrapper_config(program, ion, target, err);
+    dedx_workspace *ws;
+
     if (*err != 0)
-        return 0;
+        return 0.0f;
+    ws = allocate_wrapper_workspace(err);
+    if (*err != 0) {
+        dedx_free_config(config, &cleanup_err);
+        return 0.0f;
+    }
     dedx_load_config(ws, config, err);
+    if (*err != 0) {
+        dedx_free_config(config, &cleanup_err);
+        dedx_free_workspace(ws, &cleanup_err);
+        return 0.0f;
+    }
 
     stp = dedx_get_stp(ws, config, energy, err);
+    dedx_free_config(config, &cleanup_err);
+    dedx_free_workspace(ws, &cleanup_err);
     if (*err != 0)
-        return 0;
-    dedx_free_config(config, err);
-    dedx_free_workspace(ws, err);
+        return 0.0f;
+
     return stp;
 }
 
 int dedx_get_stp_table_size(const int program, const int ion, const int target) {
     int err = 0;
+    int cleanup_err = DEDX_OK;
     int result = -1;
-
+    dedx_config *cfg = allocate_wrapper_config(program, ion, target, &err);
     dedx_workspace *ws;
-    dedx_config *cfg = (dedx_config *) calloc(1, sizeof(dedx_config));
-    ws = dedx_allocate_workspace(1, &err);
 
-    // set program + ion + target combination
-    cfg->program = program;
-    cfg->ion = ion;
-    cfg->target = target;
+    if (err != 0)
+        return -1;
+    ws = allocate_wrapper_workspace(&err);
+    if (err != 0) {
+        dedx_free_config(cfg, &cleanup_err);
+        return -1;
+    }
 
-    // load spline database from the data files, it contains number of energy+stp samples
     dedx_load_config(ws, cfg, &err);
-
-    // config loading failed, returning exit code
-    if (err != 0)
+    if (err != 0) {
+        dedx_free_config(cfg, &cleanup_err);
+        dedx_free_workspace(ws, &cleanup_err);
         return -1;
+    }
 
-    // assume only one dataset is loaded, extract number of samples
-    if (ws->datasets == 1) {
-        result = ws->loaded_data[0]->n;
+    if (ws->active_datasets > 0) {
+        result = ws->loaded_data[cfg->cfg_id]->n;
     }
 
+    dedx_free_config(cfg, &cleanup_err);
+    dedx_free_workspace(ws, &cleanup_err);
     return result;
-};
+}
 
 int dedx_fill_default_energy_stp_table(
     const int program, const int ion, const int target, float *energies, float *stps) {
     int err = 0;
+    int cleanup_err = DEDX_OK;
     int i;
-
+    dedx_config *cfg = allocate_wrapper_config(program, ion, target, &err);
     dedx_workspace *ws;
-    dedx_config *cfg = (dedx_config *) calloc(1, sizeof(dedx_config));
-    ws = dedx_allocate_workspace(1, &err);
 
-    // set program + ion + target combination
-    cfg->program = program;
-    cfg->ion = ion;
-    cfg->target = target;
+    if (err != 0)
+        return -1;
+    ws = allocate_wrapper_workspace(&err);
+    if (err != 0) {
+        dedx_free_config(cfg, &cleanup_err);
+        return -1;
+    }
 
-    // load spline database from the data files
     dedx_load_config(ws, cfg, &err);
-
-    // config loading failed, returning exit code
-    if (err != 0)
+    if (err != 0) {
+        dedx_free_config(cfg, &cleanup_err);
+        dedx_free_workspace(ws, &cleanup_err);
         return -1;
+    }
 
-    // assume only one dataset is loaded, extract default no of energies and stopping powers
-    if (ws->datasets == 1) {
-        for (i = 0; i < ws->loaded_data[0]->n; i++) {
-            energies[i] = ws->loaded_data[0]->base[i].x;
+    if (ws->active_datasets > 0) {
+        for (i = 0; i < ws->loaded_data[cfg->cfg_id]->n; i++) {
+            energies[i] = ws->loaded_data[cfg->cfg_id]->base[i].x;
             stps[i] = dedx_get_stp(ws, cfg, energies[i], &err);
-            if (err != 0)
+            if (err != 0) {
+                dedx_free_config(cfg, &cleanup_err);
+                dedx_free_workspace(ws, &cleanup_err);
                 return -1;
+            }
         }
     }
 
+    dedx_free_config(cfg, &cleanup_err);
+    dedx_free_workspace(ws, &cleanup_err);
     return 0;
-};
+}
 
 int dedx_get_csda_range_table(const int program,
                               const int ion,
@@ -165,21 +219,34 @@ int dedx_get_csda_range_table(const int program,
                               const float *energies,
                               double *csda_ranges) {
     int err = 0;
-    dedx_config *config = (dedx_config *) calloc(1, sizeof(dedx_config));
-    config->target = target;
-    config->ion = ion;
-    config->program = program;
-    config->ion_a = _dedx_get_nucleon(config->ion, &err);
-    if (err != 0)
-        return err;
-    dedx_workspace *ws = dedx_allocate_workspace(1, &err);
+    int cleanup_err = DEDX_OK;
+    int i;
+    dedx_config *config = allocate_wrapper_config(program, ion, target, &err);
+    dedx_workspace *ws;
 
     if (err != 0)
         return err;
+    config->ion_a = dedx_internal_get_nucleon(config->ion, &err);
+    if (err != 0) {
+        dedx_free_config(config, &cleanup_err);
+        return err;
+    }
+    ws = allocate_wrapper_workspace(&err);
+    if (err != 0) {
+        dedx_free_config(config, &cleanup_err);
+        return err;
+    }
     dedx_load_config(ws, config, &err);
+    if (err != 0) {
+        dedx_free_config(config, &cleanup_err);
+        dedx_free_workspace(ws, &cleanup_err);
+        return err;
+    }
 
-    for (int i = 0; i < no_of_points; i++) {
+    for (i = 0; i < no_of_points; i++) {
         csda_ranges[i] = dedx_get_csda(ws, config, energies[i], &err);
+        if (err != 0)
+            break;
     }
     dedx_free_config(config, &err);
     dedx_free_workspace(ws, &err);
diff --git a/libdedx/dedx_wrappers.h b/libdedx/dedx_wrappers.h
index 2e4a645..d8be27c 100644
--- a/libdedx/dedx_wrappers.h
+++ b/libdedx/dedx_wrappers.h
@@ -1,5 +1,5 @@
-#ifndef DEDX_WRAPPERS_INCLUDED
-#define DEDX_WRAPPERS_INCLUDED
+#ifndef DEDX_WRAPPERS_H
+#define DEDX_WRAPPERS_H
 
 /**
  * @file dedx_wrappers.h
@@ -9,9 +9,6 @@
  * for one-off queries. For repeated evaluations use the core API in dedx.h.
  */
 
-#include <stdio.h>
-#include <stdlib.h>
-
 #include "dedx.h"
 
 /** @brief Fill an array with all supported program identifiers.
@@ -100,4 +97,4 @@ int dedx_get_csda_range_table(const int program,
                               const float *energies,
                               double *csda_ranges);
 
-#endif // DEDX_WRAPPERS_INCLUDED
+#endif // DEDX_WRAPPERS_H
diff --git a/libdedx/tools/dedx_math.h b/libdedx/tools/dedx_math.h
index e275dc7..f5730e7 100644
--- a/libdedx/tools/dedx_math.h
+++ b/libdedx/tools/dedx_math.h
@@ -1,6 +1,6 @@
-#ifndef DEDX_MATH_H_INCLUDED
-#define DEDX_MATH_H_INCLUDED
+#ifndef DEDX_MATH_H
+#define DEDX_MATH_H
 #include <stdlib.h>
 float _dedx_max(float a, float b);
 float _dedx_min(float a, float b);
-#endif
+#endif // DEDX_MATH_H
diff --git a/tests/test_core_api.c b/tests/test_core_api.c
new file mode 100644
index 0000000..e61db46
--- /dev/null
+++ b/tests/test_core_api.c
@@ -0,0 +1,145 @@
+#include <dedx_error.h>
+#include <string.h>
+
+#include "dedx_mpaul.h"
+#include "test_helpers.h"
+
+static int faili(const char *label, int got, int expected) {
+    fprintf(stderr, "FAIL %s: got %d expected %d\n", label, got, expected);
+    return 1;
+}
+
+static int failf(const char *label, double got, double expected) {
+    fprintf(stderr, "FAIL %s: got %.8g expected %.8g\n", label, got, expected);
+    return 1;
+}
+
+static int failmsg(const char *label, const char *message) {
+    fprintf(stderr, "FAIL %s: %s\n", label, message);
+    return 1;
+}
+
+static int check_string_present(const char *label, const char *value) {
+    if (value == NULL || value[0] == '\0')
+        return failmsg(label, "missing string");
+    return 0;
+}
+
+static int check_close(const char *label, double got, double expected, double rel) {
+    double scale = fabs(expected) > 1e-12 ? fabs(expected) : 1.0;
+
+    if (fabs(got - expected) > rel * scale)
+        return failf(label, got, expected);
+    return 0;
+}
+
+int main(void) {
+    int failures = 0;
+    int err = 0;
+    int major = -1;
+    int minor = -1;
+    int patch = -1;
+    const int *programs;
+    const int *materials;
+    const int *ions;
+    float composition[20][2];
+    unsigned int comp_len = 0;
+    float i_value;
+    float simple_stp;
+    dedx_workspace *ws;
+    dedx_config *cfg;
+    float direct_stp;
+    float coef_c;
+    float coef_d_low_z;
+    float coef_g;
+    float coef_h_low_z;
+    float coef_h;
+
+    dedx_get_version(&major, &minor, &patch);
+    if (major < 0 || minor < 0 || patch < 0)
+        failures += failmsg("version", "negative version component");
+
+    failures += check_string_present("program name", dedx_get_program_name(DEDX_PSTAR));
+    failures += check_string_present("program version", dedx_get_program_version(DEDX_PSTAR));
+    failures += check_string_present("material name", dedx_get_material_name(DEDX_WATER));
+    failures += check_string_present("ion name", dedx_get_ion_name(DEDX_PROTON));
+
+    programs = dedx_get_program_list();
+    if (programs == NULL || programs[0] < 0)
+        failures += failmsg("program list", "missing program list");
+
+    materials = dedx_get_material_list(DEDX_PSTAR);
+    if (materials == NULL || materials[0] < 0)
+        failures += failmsg("material list", "missing material list");
+
+    ions = dedx_get_ion_list(DEDX_DEFAULT);
+    if (ions == NULL || ions[0] != 1 || ions[111] != 112 || ions[112] != -1)
+        failures += failmsg("default ion list", "unexpected ion list contents");
+
+    dedx_get_composition(DEDX_WATER, composition, &comp_len, &err);
+    if (err != DEDX_OK) {
+        failures += faili("water composition err", err, DEDX_OK);
+    } else if (comp_len == 0) {
+        failures += failmsg("water composition", "empty composition");
+    }
+
+    err = 0;
+    i_value = dedx_get_i_value(DEDX_WATER, &err);
+    if (err != DEDX_OK) {
+        failures += faili("water i-value err", err, DEDX_OK);
+    } else if (i_value <= 0.0f) {
+        failures += failmsg("water i-value", "non-positive I value");
+    }
+
+    err = 0;
+    simple_stp = dedx_get_simple_stp(DEDX_PROTON, DEDX_WATER, 100.0f, &err);
+    if (err != DEDX_OK) {
+        failures += faili("simple stp err", err, DEDX_OK);
+    } else if (simple_stp <= 0.0f) {
+        failures += failmsg("simple stp", "non-positive stopping power");
+    }
+
+    ws = dedx_allocate_workspace(1, &err);
+    cfg = calloc(1, sizeof(dedx_config));
+    if (ws == NULL || cfg == NULL || err != DEDX_OK) {
+        fprintf(stderr, "FAIL setup direct load: err=%d\n", err);
+        dedx_free_config(cfg, &err);
+        if (ws != NULL)
+            dedx_free_workspace(ws, &err);
+        return failures + 1;
+    }
+
+    cfg->program = DEDX_ICRU;
+    cfg->ion = DEDX_PROTON;
+    cfg->target = DEDX_WATER;
+    dedx_load_config(ws, cfg, &err);
+    if (err != DEDX_OK) {
+        failures += faili("direct load err", err, DEDX_OK);
+    } else {
+        direct_stp = dedx_get_stp(ws, cfg, 100.0f, &err);
+        if (err != DEDX_OK) {
+            failures += faili("direct stp err", err, DEDX_OK);
+        } else {
+            failures += check_close("simple/direct stp", simple_stp, direct_stp, 1e-5);
+        }
+    }
+    dedx_free_config(cfg, &err);
+    dedx_free_workspace(ws, &err);
+
+    coef_c = dedx_internal_calculate_mspaul_coef('c', DEDX_CARBON, DEDX_HYDROGEN, 10.0f);
+    coef_d_low_z = dedx_internal_calculate_mspaul_coef('d', DEDX_CARBON, DEDX_HYDROGEN, 10.0f);
+    failures += check_close("mspaul d->c low-z", coef_d_low_z, coef_c, 1e-6);
+
+    coef_g = dedx_internal_calculate_mspaul_coef('g', DEDX_CARBON, DEDX_HYDROGEN, 10.0f);
+    coef_h_low_z = dedx_internal_calculate_mspaul_coef('h', DEDX_CARBON, DEDX_HYDROGEN, 10.0f);
+    failures += check_close("mspaul h->g low-z", coef_h_low_z, coef_g, 1e-6);
+
+    coef_h = dedx_internal_calculate_mspaul_coef('h', DEDX_CARBON, DEDX_WATER, 10.0f);
+    if (coef_h <= 0.0f)
+        failures += failmsg("mspaul h water", "non-positive coefficient");
+
+    if (dedx_internal_calculate_mspaul_coef('c', DEDX_HELIUM, DEDX_WATER, 10.0f) != 1.0f)
+        failures += failmsg("mspaul helium", "helium fast-path is not unity");
+
+    return failures;
+}
diff --git a/tests/test_error_codes.c b/tests/test_error_codes.c
index d9d8789..48fc96b 100644
--- a/tests/test_error_codes.c
+++ b/tests/test_error_codes.c
@@ -64,6 +64,18 @@ int main(void) {
     err = 0;
     dedx_load_config(ws, cfg, &err);
     failures += check_err(err, DEDX_ERR_RHO_REQUIRED, "BETHE custom target without rho");
+    failures += check_err(cfg->loaded, 0, "failed load should not mark config loaded");
+    failures += check_err(cfg->cfg_id, -1, "failed load should not assign cfg_id");
+    dedx_free_config(cfg, &err);
+    dedx_free_workspace(ws, &err);
+
+    /* dedx_get_stp must reject invalid dataset ids before dereferencing */
+    ws = dedx_allocate_workspace(1, &err);
+    cfg = calloc(1, sizeof(dedx_config));
+    cfg->cfg_id = 0;
+    err = 0;
+    dedx_get_stp(ws, cfg, 1.0f, &err);
+    failures += check_err(err, DEDX_ERR_INVALID_DATASET_ID, "stp with invalid dataset id");
     dedx_free_config(cfg, &err);
     dedx_free_workspace(ws, &err);
 
diff --git a/tests/test_mstar.c b/tests/test_mstar.c
index 772eef6..7781226 100644
--- a/tests/test_mstar.c
+++ b/tests/test_mstar.c
@@ -9,10 +9,119 @@ static dedx_config *make_mstar_mode_config(int target, char mode) {
     return cfg;
 }
 
+static int report_mode_equivalence_error(const char *stage, const char *label, int err) {
+    fprintf(stderr, "FAIL %s: %s err=%d\n", stage, label, err);
+    return 1;
+}
+
+static int check_mode_equivalence(int target, char lhs_mode, char rhs_mode, float energy, const char *label) {
+    int err = 0;
+    int failures = 0;
+    dedx_workspace *ws = dedx_allocate_workspace(2, &err);
+    dedx_config *lhs = make_mstar_mode_config(target, lhs_mode);
+    dedx_config *rhs = make_mstar_mode_config(target, rhs_mode);
+    float lhs_value = 0.0f;
+    float rhs_value = 0.0f;
+
+    if (err != 0 || ws == NULL || lhs == NULL || rhs == NULL) {
+        failures = report_mode_equivalence_error("alloc", label, err);
+        dedx_free_config(lhs, &err);
+        dedx_free_config(rhs, &err);
+        dedx_free_workspace(ws, &err);
+        return failures;
+    }
+
+    dedx_load_config(ws, lhs, &err);
+    if (err != 0)
+        failures = report_mode_equivalence_error("load lhs", label, err);
+
+    if (failures == 0) {
+        dedx_load_config(ws, rhs, &err);
+        if (err != 0)
+            failures = report_mode_equivalence_error("load rhs", label, err);
+    }
+
+    if (failures == 0) {
+        lhs_value = (float) dedx_get_stp(ws, lhs, energy, &err);
+        if (err != 0)
+            failures = report_mode_equivalence_error("stp lhs", label, err);
+    }
+
+    if (failures == 0) {
+        rhs_value = (float) dedx_get_stp(ws, rhs, energy, &err);
+        if (err != 0)
+            failures = report_mode_equivalence_error("stp rhs", label, err);
+    }
+
+    if (failures == 0 && check_result(lhs_value, rhs_value)) {
+        fprintf(stderr,
+                "FAIL mode-equivalence: %s target=%d E=%.3e MeV/u got %.5e expected %.5e\n",
+                label,
+                target,
+                energy,
+                lhs_value,
+                rhs_value);
+        failures = 1;
+        failures = 1;
+    }
+
+    dedx_free_config(lhs, &err);
+    dedx_free_config(rhs, &err);
+    dedx_free_workspace(ws, &err);
+    return failures;
+}
+
 int main(void) {
     int failures = 0;
     const float energy_grid[] = {0.07f, 1.0f, 10.0f, 78.0f, 1000.0f};
 
+    /* Reference values below were extracted from the original MSTAR 3.12
+     * Fortran sources via MSTAR1/MSPAUL.
+     */
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_C), energy_grid[0], 5.634276e3f, "mstar-c");
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_C), energy_grid[1], 6.592632e3f, "mstar-c");
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_C), energy_grid[2], 1.639345e3f, "mstar-c");
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_C), energy_grid[3], 3.165557e2f, "mstar-c");
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_C), energy_grid[4], 7.993779e1f, "mstar-c");
+
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_D), energy_grid[0], 5.589206e3f, "mstar-d");
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_D), energy_grid[1], 6.586625e3f, "mstar-d");
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_D), energy_grid[2], 1.639723e3f, "mstar-d");
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_D), energy_grid[3], 3.167525e2f, "mstar-d");
+    failures +=
+        check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_D), energy_grid[4], 8.000594e1f, "mstar-d");
+
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_G), energy_grid[0], 4.340127e3f, "mstar-g");
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_G), energy_grid[1], 5.262383e3f, "mstar-g");
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_G), energy_grid[2], 1.442963e3f, "mstar-g");
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_G), energy_grid[3], 2.808873e2f, "mstar-g");
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_G), energy_grid[4], 7.126015e1f, "mstar-g");
+
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_H), energy_grid[0], 4.140607e3f, "mstar-h");
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_H), energy_grid[1], 5.216270e3f, "mstar-h");
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_H), energy_grid[2], 1.447108e3f, "mstar-h");
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_H), energy_grid[3], 2.808873e2f, "mstar-h");
+    failures += check_config_stp(
+        make_mstar_mode_config(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_H), energy_grid[4], 7.126015e1f, "mstar-h");
+
     failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_WATER, energy_grid[0], 5.589e3f);
     failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_WATER, energy_grid[1], 6.587e3f);
     failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_WATER, energy_grid[2], 1.640e3f);
@@ -25,12 +134,6 @@ int main(void) {
     failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_PMMA, energy_grid[3], 3.094e2f);
     failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_PMMA, energy_grid[4], 7.762e1f);
 
-    failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_ALANINE, energy_grid[0], 6.312e3f);
-    failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_ALANINE, energy_grid[1], 6.533e3f);
-    failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_ALANINE, energy_grid[2], 1.614e3f);
-    failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_ALANINE, energy_grid[3], 3.105e2f);
-    failures += check_stp(DEDX_MSTAR, DEDX_CARBON, DEDX_ALANINE, energy_grid[4], 7.772e1f);
-
     failures +=
         check_config_stp(make_mstar_mode_config(DEDX_WATER, DEDX_MSTAR_MODE_A), energy_grid[0], 5.634e3f, "mstar-a");
     failures +=
@@ -64,5 +167,10 @@ int main(void) {
     failures +=
         check_config_stp(make_mstar_mode_config(DEDX_ALANINE, DEDX_MSTAR_MODE_C), energy_grid[4], 7.767e1f, "mstar-c");
 
+    failures +=
+        check_mode_equivalence(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_A, DEDX_MSTAR_MODE_G, 10.0f, "mstar-a-gas");
+    failures +=
+        check_mode_equivalence(DEDX_AIR_DRY_NEAR_SEA_LEVEL, DEDX_MSTAR_MODE_B, DEDX_MSTAR_MODE_H, 10.0f, "mstar-b-gas");
+
     return failures;
 }
diff --git a/tests/test_pstar.c b/tests/test_pstar.c
index 916548c..5dcf06d 100644
--- a/tests/test_pstar.c
+++ b/tests/test_pstar.c
@@ -26,7 +26,7 @@ int main(void) {
     failures += check_stp(DEDX_PSTAR, DEDX_PROTON, DEDX_ALANINE, 78.0f, 8.617e0f);
     failures += check_stp(DEDX_PSTAR, DEDX_PROTON, DEDX_ALANINE, 1000.0f, 2.148e0f);
 
-    // Proton stopping power in nitrogen (gas target — exercises _dedx_target_is_gas)
+    // Proton stopping power in nitrogen (gas target — exercises dedx_internal_target_is_gas)
     // Reference values from NIST PSTAR database
     failures += check_stp(DEDX_PSTAR, DEDX_PROTON, DEDX_NITROGEN, 0.07f, 7.624e2f);
     failures += check_stp(DEDX_PSTAR, DEDX_PROTON, DEDX_NITROGEN, 1.0f, 2.260e2f);
diff --git a/tests/test_tools.c b/tests/test_tools.c
new file mode 100644
index 0000000..b3ebc96
--- /dev/null
+++ b/tests/test_tools.c
@@ -0,0 +1,99 @@
+#include <dedx.h>
+#include <dedx_tools.h>
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+
+static int failf(const char *label, double got, double expected) {
+    fprintf(stderr, "FAIL %s: got %.8g expected %.8g\n", label, got, expected);
+    return 1;
+}
+
+static int faili(const char *label, int got, int expected) {
+    fprintf(stderr, "FAIL %s: got %d expected %d\n", label, got, expected);
+    return 1;
+}
+
+static int approx(double a, double b, double rel) {
+    double scale = fabs(b) > 1e-12 ? fabs(b) : 1.0;
+    return fabs(a - b) <= rel * scale;
+}
+
+int main(void) {
+    int failures = 0;
+    int err = 0;
+    dedx_workspace *ws = dedx_allocate_workspace(1, &err);
+    dedx_config *cfg = calloc(1, sizeof(dedx_config));
+    const float energy = 100.0f;
+    double csda;
+    double inverse_csda;
+    double stp;
+    float old_values[2] = {1.0f, 2.0f};
+    float new_values[2] = {0.0f, 0.0f};
+
+    if (ws == NULL || cfg == NULL || err != DEDX_OK) {
+        fprintf(stderr, "setup failed: err=%d\n", err);
+        dedx_free_config(cfg, &err);
+        if (ws != NULL)
+            dedx_free_workspace(ws, &err);
+        return 1;
+    }
+
+    cfg->program = DEDX_PSTAR;
+    cfg->ion = DEDX_PROTON;
+    cfg->target = DEDX_WATER;
+    cfg->ion_a = 1;
+
+    dedx_load_config(ws, cfg, &err);
+    if (err != DEDX_OK) {
+        fprintf(stderr, "dedx_load_config failed: err=%d\n", err);
+        dedx_free_config(cfg, &err);
+        dedx_free_workspace(ws, &err);
+        return 1;
+    }
+
+    csda = dedx_get_csda(ws, cfg, energy, &err);
+    if (err != DEDX_OK) {
+        failures += faili("dedx_get_csda err", err, DEDX_OK);
+    }
+
+    inverse_csda = dedx_get_inverse_csda(ws, cfg, (float) csda, &err);
+    if (err != DEDX_OK) {
+        failures += faili("dedx_get_inverse_csda err", err, DEDX_OK);
+    } else if (!approx(inverse_csda, energy, 1e-3)) {
+        failures += failf("inverse csda energy", inverse_csda, energy);
+    }
+
+    stp = dedx_get_stp(ws, cfg, energy, &err);
+    if (err != DEDX_OK) {
+        failures += faili("dedx_get_stp err", err, DEDX_OK);
+    }
+
+    err = convert_units(DEDX_MEVCM2G, DEDX_KEVUM, DEDX_WATER, 2, old_values, new_values);
+    if (err != DEDX_OK) {
+        failures += faili("convert_units err", err, DEDX_OK);
+    } else {
+        if (!approx(new_values[0], 0.1, 1e-6))
+            failures += failf("convert_units first", new_values[0], 0.1);
+        if (!approx(new_values[1], 0.2, 1e-6))
+            failures += failf("convert_units second", new_values[1], 0.2);
+    }
+
+    old_values[0] = 1.0f;
+    old_values[1] = 2.0f;
+    new_values[0] = 0.0f;
+    new_values[1] = 0.0f;
+    err = convert_units(DEDX_KEVUM, DEDX_MEVCM2G, DEDX_WATER, 2, old_values, new_values);
+    if (err != DEDX_OK) {
+        failures += faili("convert_units reverse err", err, DEDX_OK);
+    } else {
+        if (!approx(new_values[0], 10.0, 1e-6))
+            failures += failf("convert_units reverse first", new_values[0], 10.0);
+        if (!approx(new_values[1], 20.0, 1e-6))
+            failures += failf("convert_units reverse second", new_values[1], 20.0);
+    }
+
+    dedx_free_config(cfg, &err);
+    dedx_free_workspace(ws, &err);
+    return failures;
+}
diff --git a/tests/test_validate_internal.c b/tests/test_validate_internal.c
new file mode 100644
index 0000000..486a5fb
--- /dev/null
+++ b/tests/test_validate_internal.c
@@ -0,0 +1,111 @@
+#include <dedx.h>
+#include <stdio.h>
+#include <stdlib.h>
+
+#include "dedx_error.h"
+#include "dedx_validate.h"
+
+static int check_int(int got, int expected, const char *label) {
+    if (got != expected) {
+        fprintf(stderr, "FAIL %s: got %d expected %d\n", label, got, expected);
+        return 1;
+    }
+    return 0;
+}
+
+static int check_true(int condition, const char *label) {
+    if (!condition) {
+        fprintf(stderr, "FAIL %s\n", label);
+        return 1;
+    }
+    return 0;
+}
+
+static dedx_config *alloc_config(void) {
+    return calloc(1, sizeof(dedx_config));
+}
+
+static int test_validate_state_gas(void) {
+    int failures = 0;
+    int err = 0;
+    dedx_config *cfg = alloc_config();
+
+    cfg->target = DEDX_AIR;
+    cfg->compound_state = DEDX_DEFAULT_STATE;
+
+    dedx_internal_validate_state(cfg, &err);
+    failures += check_int(err, DEDX_OK, "validate_state gas err");
+    failures += check_int(cfg->compound_state, DEDX_GAS, "validate_state gas branch");
+
+    dedx_free_config(cfg, &err);
+    return failures;
+}
+
+static int test_validate_config_custom_atoms(void) {
+    int failures = 0;
+    int err = 0;
+    dedx_config *cfg = alloc_config();
+
+    cfg->program = DEDX_PSTAR;
+    cfg->target = 0;
+    cfg->elements_length = 2;
+    cfg->elements_id = calloc(2, sizeof(int));
+    cfg->elements_atoms = calloc(2, sizeof(int));
+
+    cfg->elements_id[0] = DEDX_HYDROGEN;
+    cfg->elements_id[1] = DEDX_OXYGEN;
+    cfg->elements_atoms[0] = 2;
+    cfg->elements_atoms[1] = 1;
+
+    dedx_internal_validate_config(cfg, &err);
+    failures += check_int(err, DEDX_OK, "validate_config custom atoms err");
+    failures += check_true(cfg->bragg_used == 1, "validate_config should enable Bragg rule");
+    failures += check_true(cfg->elements_mass_fraction != NULL, "validate_config should derive mass fractions");
+    if (cfg->elements_mass_fraction != NULL) {
+        failures += check_true(cfg->elements_mass_fraction[0] > 0.0f, "derived hydrogen mass fraction");
+        failures += check_true(cfg->elements_mass_fraction[1] > cfg->elements_mass_fraction[0],
+                               "oxygen mass fraction should dominate water-like mixture");
+    }
+
+    dedx_free_config(cfg, &err);
+    return failures;
+}
+
+static int test_evaluate_i_pot_custom_elements(void) {
+    int failures = 0;
+    int err = 0;
+    dedx_config *cfg = alloc_config();
+
+    cfg->target = 0;
+    cfg->compound_state = DEDX_CONDENSED;
+    cfg->elements_length = 2;
+    cfg->elements_id = calloc(2, sizeof(int));
+    cfg->elements_mass_fraction = calloc(2, sizeof(float));
+
+    cfg->elements_id[0] = DEDX_HYDROGEN;
+    cfg->elements_id[1] = DEDX_OXYGEN;
+    cfg->elements_mass_fraction[0] = 0.111894f;
+    cfg->elements_mass_fraction[1] = 0.888106f;
+
+    dedx_internal_evaluate_i_pot(cfg, &err);
+    failures += check_int(err, DEDX_OK, "evaluate_i_pot custom err");
+    failures += check_true(cfg->elements_i_value != NULL, "evaluate_i_pot should allocate element I values");
+    if (cfg->elements_i_value != NULL) {
+        failures += check_true(cfg->elements_i_value[0] > 0.0f, "hydrogen I value");
+        failures += check_true(cfg->elements_i_value[1] > 0.0f, "oxygen I value");
+    }
+    failures += check_true(cfg->i_value > 0.0f, "compound I value");
+
+    dedx_free_config(cfg, &err);
+    return failures;
+}
+
+int main(void) {
+    int failures = 0;
+
+    failures += test_validate_state_gas();
+    failures += test_validate_config_custom_atoms();
+    failures += test_evaluate_i_pot_custom_elements();
+
+    return failures;
+}