Move sources into src

7 files changed, 1045 insertions, 0 deletions

diff --git a/src/deco.c b/src/deco.c
new file mode 100644
index 0000000..5ef6537
--- /dev/null
+++ b/src/deco.c
@@ -0,0 +1,306 @@
+/* SPDX-License-Identifier: MIT-0 */
+
+#include <assert.h>
+#include <math.h>
+#include <stdbool.h>
+
+#include "deco.h"
+
+enum ALGO ALGO_VER = ZHL_16C;
+
+#define PO2_MAX (1.6)
+#define END_MAX (abs_depth(msw_to_bar(30)))
+#define RND(x) (round((x) *10000) / 10000)
+
+typedef struct zhl_n2_t {
+    double t;
+    double a[3];
+    double b;
+} zhl_n2_t;
+
+typedef struct zhl_he_t {
+    double t;
+    double a;
+    double b;
+} zhl_he_t;
+
+const zhl_n2_t ZHL16N[] = {
+    {.t = 5.0,   .a = {1.1696, 1.1696, 1.1696}, .b = 0.5578},
+    {.t = 8.0,   .a = {1.0000, 1.0000, 1.0000}, .b = 0.6514},
+    {.t = 12.5,  .a = {0.8618, 0.8618, 0.8618}, .b = 0.7222},
+    {.t = 18.5,  .a = {0.7562, 0.7562, 0.7562}, .b = 0.7825},
+    {.t = 27.0,  .a = {0.6667, 0.6667, 0.6200}, .b = 0.8126},
+    {.t = 38.3,  .a = {0.5933, 0.5600, 0.5043}, .b = 0.8434},
+    {.t = 54.3,  .a = {0.5282, 0.4947, 0.4410}, .b = 0.8693},
+    {.t = 77.0,  .a = {0.4701, 0.4500, 0.4000}, .b = 0.8910},
+    {.t = 109.0, .a = {0.4187, 0.4187, 0.3750}, .b = 0.9092},
+    {.t = 146.0, .a = {0.3798, 0.3798, 0.3500}, .b = 0.9222},
+    {.t = 187.0, .a = {0.3497, 0.3497, 0.3295}, .b = 0.9319},
+    {.t = 239.0, .a = {0.3223, 0.3223, 0.3065}, .b = 0.9403},
+    {.t = 305.0, .a = {0.2971, 0.2850, 0.2835}, .b = 0.9477},
+    {.t = 390.0, .a = {0.2737, 0.2737, 0.2610}, .b = 0.9544},
+    {.t = 498.0, .a = {0.2523, 0.2523, 0.2480}, .b = 0.9602},
+    {.t = 635.0, .a = {0.2327, 0.2327, 0.2327}, .b = 0.9653},
+};
+
+const zhl_he_t ZHL16He[] = {
+    {.t = 1.88,   .a = 1.6189, .b = 0.4770},
+    {.t = 3.02,   .a = 1.3830, .b = 0.5747},
+    {.t = 4.72,   .a = 1.1919, .b = 0.6527},
+    {.t = 6.99,   .a = 1.0458, .b = 0.7223},
+    {.t = 10.21,  .a = 0.9220, .b = 0.7582},
+    {.t = 14.48,  .a = 0.8205, .b = 0.7957},
+    {.t = 20.53,  .a = 0.7305, .b = 0.8279},
+    {.t = 29.11,  .a = 0.6502, .b = 0.8553},
+    {.t = 41.20,  .a = 0.5950, .b = 0.8757},
+    {.t = 55.19,  .a = 0.5545, .b = 0.8903},
+    {.t = 70.69,  .a = 0.5333, .b = 0.8997},
+    {.t = 90.34,  .a = 0.5189, .b = 0.9073},
+    {.t = 115.29, .a = 0.5181, .b = 0.9122},
+    {.t = 147.42, .a = 0.5176, .b = 0.9171},
+    {.t = 188.24, .a = 0.5172, .b = 0.9217},
+    {.t = 240.03, .a = 0.5119, .b = 0.9267},
+};
+
+double bar_to_msw(const double bar)
+{
+    return bar * 10;
+}
+
+double msw_to_bar(const double msw)
+{
+    return msw / 10;
+}
+
+double abs_depth(const double gd)
+{
+    return gd + SURFACE_PRESSURE;
+}
+
+double gauge_depth(const double ad)
+{
+    return ad - SURFACE_PRESSURE;
+}
+
+gas_t gas_new(const unsigned char o2, const unsigned char he, double mod)
+{
+    assert(o2 + he <= 100);
+
+    if (mod == MOD_AUTO) {
+        double mod_po2 = PO2_MAX / (o2 / 100.0);
+        double mod_end = END_MAX / (1 - he / 100.0);
+
+        mod = min(mod_po2, mod_end);
+    }
+
+    return (gas_t){.o2 = o2, .he = he, .n2 = 100 - o2 - he, .mod = mod};
+}
+
+int gas_equal(const gas_t *g1, const gas_t *g2)
+{
+    return g1->o2 == g2->o2 && g1->he == g2->he && g1->mod == g2->mod;
+}
+
+unsigned char gas_o2(const gas_t *gas)
+{
+    return gas->o2;
+}
+
+unsigned char gas_he(const gas_t *gas)
+{
+    return gas->he;
+}
+
+unsigned char gas_n2(const gas_t *gas)
+{
+    return gas->n2;
+}
+
+double gas_mod(const gas_t *gas)
+{
+    return gas->mod;
+}
+
+double add_segment_ascdec(decostate_t *ds, const double dstart, const double dend, const double time, const gas_t *gas)
+{
+    assert(time > 0);
+
+    const double rate = (dend - dstart) / time;
+
+    for (int i = 0; i < 16; i++) {
+        double pio = gas_he(gas) / 100.0 * (dstart - P_WV);
+        double po = ds->phe[i];
+        double r = gas_he(gas) / 100.0 * rate;
+        double k = log(2) / ZHL16He[i].t;
+        double t = time;
+
+        ds->phe[i] = pio + r * (t - 1 / k) - (pio - po - (r / k)) * exp(-k * t);
+    }
+
+    for (int i = 0; i < 16; i++) {
+        double pio = gas_n2(gas) / 100.0 * (dstart - P_WV);
+        double po = ds->pn2[i];
+        double r = gas_n2(gas) / 100.0 * rate;
+        double k = log(2) / ZHL16N[i].t;
+        double t = time;
+
+        ds->pn2[i] = pio + r * (t - 1 / k) - (pio - po - (r / k)) * exp(-k * t);
+    }
+
+    /* TODO add CNS */
+    /* TODO add OTU */
+
+    if (dend > ds->max_depth)
+        ds->max_depth = dend;
+
+    return time;
+}
+
+double add_segment_const(decostate_t *ds, const double depth, const double time, const gas_t *gas)
+{
+    assert(time > 0);
+
+    for (int i = 0; i < 16; i++) {
+        double pio = gas_he(gas) / 100.0 * (depth - P_WV);
+        double po = ds->phe[i];
+        double k = log(2) / ZHL16He[i].t;
+        double t = time;
+
+        ds->phe[i] = po + (pio - po) * (1 - exp(-k * t));
+    }
+
+    for (int i = 0; i < 16; i++) {
+        double pio = gas_n2(gas) / 100.0 * (depth - P_WV);
+        double po = ds->pn2[i];
+        double k = log(2) / ZHL16N[i].t;
+        double t = time;
+
+        ds->pn2[i] = po + (pio - po) * (1 - exp(-k * t));
+    }
+
+    /* TODO add CNS */
+    /* TODO add OTU */
+
+    if (depth > ds->max_depth)
+        ds->max_depth = depth;
+
+    return time;
+}
+
+double get_gf(const decostate_t *ds, const double depth)
+{
+    const unsigned char lo = ds->gflo;
+    const unsigned char hi = ds->gfhi;
+
+    if (ds->firststop == -1)
+        return lo;
+
+    if (depth < SURFACE_PRESSURE)
+        return hi;
+
+    if (depth > ds->firststop)
+        return lo;
+
+    /* interpolate lo and hi between first stop and surface */
+    double next_stop = depth - ds->ceil_multiple;
+    return hi - (hi - lo) * gauge_depth(next_stop) / gauge_depth(ds->firststop);
+}
+
+double round_ceiling(const decostate_t *ds, const double c)
+{
+    assert(ds->ceil_multiple != 0);
+
+    return abs_depth(ds->ceil_multiple * ceil(RND(gauge_depth(c) / ds->ceil_multiple)));
+}
+
+double ceiling(const decostate_t *ds, double gf)
+{
+    double c = 0;
+    gf /= 100;
+
+    for (int i = 0; i < 16; i++) {
+        /* n2 a and b values */
+        double an = ZHL16N[i].a[ALGO_VER];
+        double bn = ZHL16N[i].b;
+
+        /* he a and b values */
+        double ah = ZHL16He[i].a;
+        double bh = ZHL16He[i].b;
+
+        /* scale n2 and he values for a and b proportional to their pressure */
+        double pn2 = ds->pn2[i];
+        double phe = ds->phe[i];
+
+        double a = ((an * pn2) + (ah * phe)) / (pn2 + phe);
+        double b = ((bn * pn2) + (bh * phe)) / (pn2 + phe);
+
+        /* update ceiling */
+        c = max(c, ((pn2 + phe) - (a * gf)) / (gf / b + 1 - gf));
+    }
+
+    return round_ceiling(ds, c);
+}
+
+double gf99(const decostate_t *ds, double depth)
+{
+    double gf = 0;
+
+    for (int i = 0; i < 16; i++) {
+        /* n2 a and b values */
+        double an = ZHL16N[i].a[ALGO_VER];
+        double bn = ZHL16N[i].b;
+
+        /* he a and b values */
+        double ah = ZHL16He[i].a;
+        double bh = ZHL16He[i].b;
+
+        /* scale n2 and he values for a and b proportional to their pressure */
+        double pn2 = ds->pn2[i];
+        double phe = ds->phe[i];
+
+        double a = ((an * pn2) + (ah * phe)) / (pn2 + phe);
+        double b = ((bn * pn2) + (bh * phe)) / (pn2 + phe);
+
+        /* update gf99 */
+        gf = max(gf, (pn2 + phe - depth) / (a + depth / b - depth));
+    }
+
+    return gf * 100;
+}
+
+void init_tissues(decostate_t *ds)
+{
+    const double pn2 = 0.79 * (SURFACE_PRESSURE - P_WV);
+    const double phe = 0.00 * (SURFACE_PRESSURE - P_WV);
+
+    for (int i = 0; i < 16; i++)
+        ds->pn2[i] = pn2;
+
+    for (int i = 0; i < 16; i++)
+        ds->phe[i] = phe;
+}
+
+void init_decostate(decostate_t *ds, const unsigned char gflo, const unsigned char gfhi, const double ceil_multiple)
+{
+    init_tissues(ds);
+
+    ds->gflo = gflo;
+    ds->gfhi = gfhi;
+    ds->firststop = -1;
+    ds->ceil_multiple = ceil_multiple;
+}
+
+double ppO2(double depth, const gas_t *gas)
+{
+    return gas_o2(gas) / 100.0 * depth;
+}
+
+double end(double depth, const gas_t *gas)
+{
+    return (gas_o2(gas) + gas_n2(gas)) / 100.0 * depth;
+}
+
+double ead(double depth, const gas_t *gas)
+{
+    return depth * gas_n2(gas) / 79.0;
+}
diff --git a/src/deco.h b/src/deco.h
new file mode 100644
index 0000000..439626f
--- /dev/null
+++ b/src/deco.h
@@ -0,0 +1,69 @@
+/* SPDX-License-Identifier: MIT-0 */
+
+#ifndef DECO_H
+#define DECO_H
+
+#include <stddef.h>
+
+#define max(X, Y) (((X) > (Y)) ? (X) : (Y))
+#define min(X, Y) (((X) < (Y)) ? (X) : (Y))
+#define len(X) (sizeof(X) / sizeof((X)[0]))
+
+#define P_WV_BUHL 0.0627 /* Buhlmann value, Rq = 1.0, least conservative */
+#define P_WV_NAVY 0.0567 /* US. Navy value, Rq = 0.9 */
+#define P_WV_SCHR 0.0493 /* Schreiner value, Rq = 0.8, most conservative */
+#define P_WV P_WV_BUHL
+
+#define SURFACE_PRESSURE 1.01325
+#define MOD_AUTO 0
+
+enum ALGO {
+    ZHL_16A = 0,
+    ZHL_16B = 1,
+    ZHL_16C = 2,
+};
+
+typedef struct decostate_t {
+    double pn2[16];
+    double phe[16];
+    unsigned char gflo;
+    unsigned char gfhi;
+    double firststop;
+    double max_depth;
+    double ceil_multiple;
+} decostate_t;
+
+typedef struct gas_t {
+    unsigned char o2;
+    unsigned char he;
+    unsigned char n2;
+    double mod;
+} gas_t;
+
+double bar_to_msw(const double bar);
+double msw_to_bar(const double msw);
+double abs_depth(const double gd);
+double gauge_depth(const double ad);
+
+gas_t gas_new(const unsigned char o2, const unsigned char he, double mod);
+int gas_equal(const gas_t *g1, const gas_t *g2);
+unsigned char gas_o2(const gas_t *gas);
+unsigned char gas_he(const gas_t *gas);
+unsigned char gas_n2(const gas_t *gas);
+double gas_mod(const gas_t *gas);
+
+double add_segment_ascdec(decostate_t *ds, const double dstart, const double dend, const double time,
+                          const gas_t *gas);
+double add_segment_const(decostate_t *ds, const double depth, const double time, const gas_t *gas);
+double get_gf(const decostate_t *ds, const double depth);
+double round_ceiling(const decostate_t *ds, const double c);
+double ceiling(const decostate_t *ds, double gf);
+double gf99(const decostate_t *ds, double depth);
+
+void init_decostate(decostate_t *ds, const unsigned char gflo, const unsigned char gfhi, const double ceil_multiple);
+
+double ppO2(double depth, const gas_t *gas);
+double end(double depth, const gas_t *gas);
+double ead(double depth, const gas_t *gas);
+
+#endif /* end of include guard: DECO_H */
diff --git a/src/opendeco.c b/src/opendeco.c
new file mode 100644
index 0000000..9231ab3
--- /dev/null
+++ b/src/opendeco.c
@@ -0,0 +1,208 @@
+/* SPDX-License-Identifier: MIT-0 */
+
+#include <argp.h>
+#include <locale.h>
+#include <math.h>
+#include <stdlib.h>
+#include <string.h>
+#include <wchar.h>
+
+#include "deco.h"
+#include "schedule.h"
+#include "output.h"
+
+#define MOD_OXY (abs_depth(msw_to_bar(6)))
+
+#ifndef VERSION
+#define VERSION "unknown version"
+#endif
+
+/* argp settings */
+static char args_doc[] = "";
+static char doc[] = "Implementation of Buhlmann ZH-L16 with Gradient Factors:"
+                    "\vExamples:\n\n"
+                    "\t./opendeco -d 18 -t 60 -g Air\n"
+                    "\t./opendeco -d 30 -t 60 -g EAN32\n"
+                    "\t./opendeco -d 40 -t 120 -g 21/35 -l 20 -h 80 --decogasses Oxygen,EAN50\n";
+const char *argp_program_bug_address = "<~tsegers/opendeco@lists.sr.ht> or https://todo.sr.ht/~tsegers/opendeco";
+const char *argp_program_version = "opendeco " VERSION;
+
+static struct argp_option options[] = {
+    {"depth",      'd', "NUMBER", 0, "Set the depth of the dive in meters",                       0},
+    {"time",       't', "NUMBER", 0, "Set the time of the dive in minutes",                       1},
+    {"gas",        'g', "STRING", 0, "Set the bottom gas used during the dive, defaults to Air",  2},
+    {"gflow",      'l', "NUMBER", 0, "Set the gradient factor at the first stop, defaults to 30", 3},
+    {"gfhigh",     'h', "NUMBER", 0, "Set the gradient factor at the surface, defaults to 75",    4},
+    {"decogasses", 'G', "LIST",   0, "Set the gasses available for deco",                         5},
+    {0,            0,   0,        0, 0,                                                           0}
+};
+
+struct arguments {
+    double depth;
+    double time;
+    char *gas;
+    int gflow;
+    int gfhigh;
+    char *decogasses;
+};
+
+static error_t parse_opt(int key, char *arg, struct argp_state *state)
+{
+    struct arguments *arguments = state->input;
+
+    switch (key) {
+    case 'd':
+        arguments->depth = arg ? atof(arg) : -1;
+        break;
+    case 't':
+        arguments->time = arg ? atof(arg) : -1;
+        break;
+    case 'g':
+        arguments->gas = arg;
+        break;
+    case 'G':
+        arguments->decogasses = arg;
+        break;
+    case 'l':
+        arguments->gflow = arg ? atoi(arg) : 100;
+        break;
+    case 'h':
+        arguments->gfhigh = arg ? atoi(arg) : 100;
+        break;
+    case ARGP_KEY_END:
+        if (arguments->depth < 0 || arguments->time < 0) {
+            argp_state_help(state, stderr, ARGP_HELP_USAGE);
+            argp_failure(state, 1, 0, "Options -d and -t are required. See --help for more information");
+            exit(ARGP_ERR_UNKNOWN);
+        }
+    default:
+        return ARGP_ERR_UNKNOWN;
+    }
+
+    return 0;
+}
+
+static struct argp argp = {options, parse_opt, args_doc, doc, 0, 0, 0};
+
+void print_segment_callback(const decostate_t *ds, const waypoint_t wp, segtype_t type)
+{
+    static double last_depth;
+    static double runtime;
+
+    wchar_t sign;
+
+    runtime += wp.time;
+
+    if (wp.depth < last_depth)
+        sign = ASC;
+    else if (wp.depth > last_depth)
+        sign = DEC;
+    else
+        sign = LVL;
+
+    if (type != SEG_TRAVEL)
+        print_planline(sign, wp.depth, wp.time, runtime, wp.gas);
+
+    last_depth = wp.depth;
+}
+
+int parse_gasses(gas_t **gasses, char *str)
+{
+    if (!str) {
+        *gasses = NULL;
+        return 0;
+    }
+
+    /* count number of gasses in string */
+    int nof_gasses = 1;
+
+    for (int c = 0; str[c]; c++)
+        if (str[c] == ',')
+            nof_gasses++;
+
+    /* allocate gas array */
+    gas_t *deco_gasses = malloc(nof_gasses * sizeof(gas_t));
+
+    /* fill gas array */
+    char *gas_str = NULL;
+    int gas_idx = 0;
+
+    while (1) {
+        if (!gas_str)
+            gas_str = strtok(str, ",");
+        else
+            gas_str = strtok(NULL, ",");
+
+        if (!gas_str)
+            break;
+
+        scan_gas(&deco_gasses[gas_idx], gas_str);
+        gas_idx++;
+    }
+
+    *gasses = deco_gasses;
+    return nof_gasses;
+}
+
+int main(int argc, char *argv[])
+{
+    setlocale(LC_ALL, "en_US.utf8");
+
+    /* argp */
+    struct arguments arguments;
+
+    arguments.depth = -1;
+    arguments.time = -1;
+    arguments.gas = "Air";
+    arguments.gflow = 30;
+    arguments.gfhigh = 75;
+    arguments.decogasses = "";
+
+    argp_parse(&argp, argc, argv, 0, 0, &arguments);
+
+    /* setup */
+    decostate_t ds;
+    init_decostate(&ds, arguments.gflow, arguments.gfhigh, msw_to_bar(3));
+    double dec_per_min = msw_to_bar(9);
+
+    gas_t bottom_gas;
+    scan_gas(&bottom_gas, arguments.gas);
+
+    gas_t *deco_gasses;
+    int nof_gasses = parse_gasses(&deco_gasses, arguments.decogasses);
+
+    /* override oxygen mod */
+    for (int i = 0; i < nof_gasses; i++)
+        if (gas_o2(&deco_gasses[i]) == 100)
+            deco_gasses[i].mod = MOD_OXY;
+
+    /* simulate dive */
+    double descent_time = msw_to_bar(arguments.depth) / dec_per_min;
+    double bottom_time = max(1, arguments.time - descent_time);
+
+    waypoint_t waypoints[] = {
+        {.depth = abs_depth(msw_to_bar(arguments.depth)), .time = descent_time, &bottom_gas},
+        {.depth = abs_depth(msw_to_bar(arguments.depth)), .time = bottom_time,  &bottom_gas},
+    };
+
+    print_planhead();
+    simulate_dive(&ds, waypoints, len(waypoints), &print_segment_callback);
+
+    /* generate deco schedule */
+    double depth = waypoints[len(waypoints) - 1].depth;
+    const gas_t *gas = waypoints[len(waypoints) - 1].gas;
+
+    /* determine @+5 TTS */
+    decostate_t ds_ = ds;
+    add_segment_const(&ds_, depth, 5, gas);
+    decoinfo_t di_plus5 = calc_deco(&ds_, depth, gas, deco_gasses, nof_gasses, NULL);
+
+    /* print actual deco schedule */
+    decoinfo_t di = calc_deco(&ds, depth, gas, deco_gasses, nof_gasses, &print_segment_callback);
+
+    /* output deco info and disclaimer */
+    wprintf(L"\nNDL: %i TTS: %i TTS @+5: %i\n", (int) floor(di.ndl), (int) ceil(di.tts), (int) ceil(di_plus5.tts));
+    print_planfoot(&ds);
+
+    return 0;
+}
diff --git a/src/output.c b/src/output.c
new file mode 100644
index 0000000..589c811
--- /dev/null
+++ b/src/output.c
@@ -0,0 +1,129 @@
+/* SPDX-License-Identifier: MIT-0 */
+
+#include <math.h>
+#include <stdio.h>
+#include <string.h>
+
+#include "output.h"
+
+extern enum ALGO ALGO_VER;
+
+void format_mm_ss(char *buf, const size_t buflen, const double time)
+{
+    double mm;
+    double ss = round(modf(time, &mm) * 60);
+
+    /* prevents 0.99999 minutes showing as 00:60 */
+    mm += ss / 60;
+    ss = (int) ss % 60;
+
+    snprintf(buf, buflen, "%3i:%02i", (int) mm, (int) ss);
+}
+
+void format_gas(char *buf, const size_t buflen, const gas_t *gas)
+{
+    if (gas_o2(gas) == 21 && gas_he(gas) == 0)
+        snprintf(buf, buflen, "Air");
+    else if (gas_o2(gas) == 100)
+        snprintf(buf, buflen, "Oxygen");
+    else if (gas_he(gas) == 0)
+        snprintf(buf, buflen, "Nitrox %i", gas_o2(gas));
+    else
+        snprintf(buf, buflen, "%i/%i", gas_o2(gas), gas_he(gas));
+}
+
+void scan_gas(gas_t *gas, char *str)
+{
+    int o2 = 0;
+    int he = 0;
+
+    if (!strcmp(str, "Air")) {
+        *gas = gas_new(21, 0, MOD_AUTO);
+        return;
+    } else if (!strcmp(str, "Oxygen")) {
+        *gas = gas_new(100, 0, MOD_AUTO);
+        return;
+    } else if (!strncmp(str, "EAN", strlen("EAN"))) {
+        sscanf(str, "EAN%i", &o2);
+    } else if (!strncmp(str, "Nitrox", strlen("Nitrox"))) {
+        sscanf(str, "Nitrox %i", &o2);
+    } else {
+        sscanf(str, "%i/%i", &o2, &he);
+    }
+
+    *gas = gas_new(o2, he, MOD_AUTO);
+}
+
+void print_planhead()
+{
+    wprintf(L"DIVE PLAN\n\n");
+    wprintf(L" %-1s  %-5s  %-8s  %-7s  %1s %-9s  %-4s  %-3s\n", "", "Depth", "Duration", "Runtime", "", "Gas", "pO2",
+            "EAD");
+}
+
+void print_planline(const wchar_t sign, const double depth, const double time, const double runtime, const gas_t *gas)
+{
+    static char gasbuf[11];
+    static char runbuf[8];
+    static char pO2buf[5];
+    static char eadbuf[4];
+    static char timbuf[16];
+
+    static gas_t last_gas;
+
+    const int depth_m = round(bar_to_msw(gauge_depth(depth)));
+    const int ead_m = round(bar_to_msw(max(0, gauge_depth(ead(depth, gas)))));
+
+    wchar_t swi = L' ';
+
+    snprintf(runbuf, len(runbuf), "(%i)", (int) ceil(runtime));
+    format_gas(gasbuf, len(gasbuf), gas);
+    format_mm_ss(timbuf, len(timbuf), time);
+
+    /* print gas swich symbol if gas changed */
+    if (!gas_equal(gas, &last_gas)) {
+        last_gas = *gas;
+        swi = SWI;
+    }
+
+    /* only print ead and pO2 on stops */
+    if (sign == LVL) {
+        snprintf(eadbuf, 4, "%3i", ead_m);
+        snprintf(pO2buf, 5, "%4.2f", ppO2(depth, gas));
+    } else {
+        snprintf(eadbuf, 4, "%3s", "-");
+        snprintf(pO2buf, 5, "%4s", "-");
+    }
+
+    wprintf(L" %lc  %4im  %8s  %-7s  %lc %-9s  %s  %s\n", sign, depth_m, timbuf, runbuf, swi, gasbuf, pO2buf, eadbuf);
+}
+
+void print_planfoot(const decostate_t *ds)
+{
+    char *model;
+    char *rq;
+
+    if (ALGO_VER == ZHL_16A)
+        model = "ZHL-16A";
+    else if (ALGO_VER == ZHL_16B)
+        model = "ZHL-16B";
+    else if (ALGO_VER == ZHL_16C)
+        model = "ZHL-16C";
+    else
+        model = "???";
+
+    if (P_WV == P_WV_BUHL)
+        rq = "1.0";
+    else if (P_WV == P_WV_NAVY)
+        rq = "0.9";
+    else if (P_WV == P_WV_SCHR)
+        rq = "0.8";
+    else
+        rq = "???";
+
+    wprintf(L"\nDeco model: Buhlmann %s\n", model);
+    wprintf(L"Conservatism: GF %i/%i, Rq = %s\n", ds->gflo, ds->gfhi, rq);
+    wprintf(L"Surface pressure: %4.3fbar\n\n", SURFACE_PRESSURE);
+
+    wprintf(L"WARNING: DIVE PLAN MAY BE INACCURATE AND MAY CONTAIN\nERRORS THAT COULD LEAD TO INJURY OR DEATH.\n");
+}
diff --git a/src/output.h b/src/output.h
new file mode 100644
index 0000000..640b457
--- /dev/null
+++ b/src/output.h
@@ -0,0 +1,21 @@
+/* SPDX-License-Identifier: MIT-0 */
+
+#ifndef OUTPUT_H
+#define OUTPUT_H
+
+#include <wchar.h>
+
+#include "deco.h"
+
+#define ASC 0x2197 /* Unicode North East Arrow */
+#define LVL 0x2192 /* Unicode Rightwards Arrow */
+#define DEC 0x2198 /* Unicode South East Arrow */
+#define SWI 0x21BB /* Clockwise Open Circle Arrow */
+
+void print_planhead();
+void print_planline(const wchar_t sign, const double depth, const double time, const double runtime, const gas_t *gas);
+void print_planfoot(const decostate_t *ds);
+
+void scan_gas(gas_t *gas, char *str);
+
+#endif /* end of include guard: OUTPUT_H */
diff --git a/src/schedule.c b/src/schedule.c
new file mode 100644
index 0000000..d13a9c9
--- /dev/null
+++ b/src/schedule.c
@@ -0,0 +1,270 @@
+/* SPDX-License-Identifier: MIT-0 */
+
+#include <assert.h>
+#include <math.h>
+
+#include "schedule.h"
+
+#define SWITCH_INTERMEDIATE 1
+
+#define STOPLEN_ROUGH 10
+#define STOPLEN_FINE 1
+
+const gas_t *best_gas(const double depth, const gas_t *gasses, const int nof_gasses)
+{
+    const gas_t *best = NULL;
+    double mod_best = -1;
+
+    for (int i = 0; i < nof_gasses; i++) {
+        double mod = gas_mod(&gasses[i]);
+
+        if (depth <= mod && (mod_best == -1 || mod < mod_best)) {
+            best = &gasses[i];
+            mod_best = mod;
+        }
+    }
+
+    return best;
+}
+
+const gas_t *next_gas(const double depth, const gas_t *gasses, const int nof_gasses)
+{
+    const gas_t *next = NULL;
+    double mod_best = 0;
+
+    for (int i = 0; i < nof_gasses; i++) {
+        double mod = gas_mod(&gasses[i]);
+
+        if (depth > mod && mod > mod_best) {
+            next = &gasses[i];
+            mod_best = mod;
+        }
+    }
+
+    return next;
+}
+
+int direct_ascent(const decostate_t *ds, const double depth, const double time, const gas_t *gas)
+{
+    decostate_t ds_ = *ds;
+    assert(ds_.firststop == -1);
+
+    add_segment_ascdec(&ds_, depth, abs_depth(0), time, gas);
+
+    return gauge_depth(ceiling(&ds_, ds_.gfhi)) <= 0;
+}
+
+void simulate_dive(decostate_t *ds, waypoint_t *waypoints, const int nof_waypoints, waypoint_callback_t wp_cb)
+{
+    double depth = abs_depth(0);
+    double runtime = 0;
+
+    for (int i = 0; i < nof_waypoints; i++) {
+        double d = waypoints[i].depth;
+        double t = waypoints[i].time;
+        const gas_t *g = waypoints[i].gas;
+
+        if (d != depth) {
+            runtime += add_segment_ascdec(ds, depth, d, t, g);
+            depth = d;
+        } else {
+            runtime += add_segment_const(ds, d, t, g);
+        }
+
+        if (wp_cb)
+            wp_cb(ds, (waypoint_t){.depth = d, .time = t, .gas = g}, SEG_DIVE);
+    }
+}
+
+double calc_ndl(decostate_t *ds, const double depth, const double ascrate, const gas_t *gas)
+{
+    decostate_t ds_ = *ds;
+    double ndl = 0;
+
+    while (ndl < 360) {
+        double tmp = add_segment_const(&ds_, depth, 1, gas);
+
+        if (!direct_ascent(&ds_, depth, gauge_depth(depth) / ascrate, gas))
+            break;
+
+        ndl += tmp;
+    }
+
+    return ndl;
+}
+
+double calc_stoplen_rough(const decostate_t *ds, const double depth, const double current_gf, const gas_t *gas)
+{
+    decostate_t ds_ = *ds;
+    double stoplen = 0;
+
+    for (;;) {
+        double tmp = add_segment_const(&ds_, depth, STOPLEN_ROUGH, gas);
+
+        if (ceiling(&ds_, current_gf) != depth)
+            break;
+
+        stoplen += tmp;
+    }
+
+    return stoplen;
+}
+
+double deco_stop(decostate_t *ds, const double depth, const double current_gf, const gas_t *gas)
+{
+    double stoplen = 0;
+
+    /* rough steps */
+    double stoplen_rough = calc_stoplen_rough(ds, depth, current_gf, gas);
+
+    if (stoplen_rough) {
+        add_segment_const(ds, depth, stoplen_rough, gas);
+        stoplen += stoplen_rough;
+    }
+
+    /* fine steps */
+    while (ceiling(ds, current_gf) == depth)
+        stoplen += add_segment_const(ds, depth, STOPLEN_FINE, gas);
+
+    return stoplen;
+}
+
+int should_switch(const gas_t **next, double *switch_depth, const decostate_t *ds, const double depth,
+                  const double next_stop, const gas_t *deco_gasses, const int nof_gasses)
+{
+    /* check if we switch at MOD or just at stops */
+    if (!SWITCH_INTERMEDIATE)
+        return 0;
+
+    /* check if there is a gas to switch to */
+    *next = next_gas(depth, deco_gasses, nof_gasses);
+
+    if (*next == NULL)
+        return 0;
+
+    /* check if the switch happens before the current ceiling */
+    *switch_depth = round_ceiling(ds, (*next)->mod) - ds->ceil_multiple;
+    assert(*switch_depth <= (*next)->mod);
+
+    if (*switch_depth <= next_stop)
+        return 0;
+
+    return 1;
+}
+
+decoinfo_t calc_deco(decostate_t *ds, const double start_depth, const gas_t *start_gas, const gas_t *deco_gasses,
+                     const int nof_gasses, waypoint_callback_t wp_cb)
+{
+    decoinfo_t ret = {.tts = 0, .ndl = 0};
+
+    /* setup start parameters */
+    double depth = start_depth;
+    const gas_t *gas = start_gas;
+
+    const double asc_per_min = msw_to_bar(9);
+
+    /* check if direct ascent is possible */
+    if (direct_ascent(ds, depth, gauge_depth(depth) / asc_per_min, gas)) {
+        ret.ndl = calc_ndl(ds, depth, asc_per_min, gas);
+        return ret;
+    }
+
+    /* determine first stop */
+    double current_gf = get_gf(ds, depth);
+
+    if (ds->firststop == -1)
+        ds->firststop = ceiling(ds, current_gf);
+
+    /* switch to best deco gas if there is one available */
+    const gas_t *best = best_gas(depth, deco_gasses, nof_gasses);
+
+    if (best)
+        gas = best;
+
+    /* alternate between ascending and stopping */
+    for (;;) {
+        /* extra bookkeeping because waypoints and segments do not match 1:1 */
+        double last_waypoint_depth = depth;
+        double waypoint_time;
+
+        /* ascend */
+        for (;;) {
+            /* determine next stop */
+            double next_stop = ceiling(ds, current_gf);
+
+            /* find out if we need to switch gas on the way */
+            const gas_t *next;
+            double switch_depth;
+
+            if (should_switch(&next, &switch_depth, ds, depth, next_stop, deco_gasses, nof_gasses)) {
+                /* ascend to gas switch */
+                ret.tts += add_segment_ascdec(ds, depth, switch_depth, fabs(depth - switch_depth) / asc_per_min, gas);
+                depth = switch_depth;
+                current_gf = get_gf(ds, depth);
+
+                /*
+                 * since we're stopping for a switch next, this is the last of
+                 * any number of consecutive travel segments and the waypoint
+                 * callback should be called.
+                 */
+                waypoint_time = fabs(last_waypoint_depth - depth) / asc_per_min;
+
+                if (wp_cb)
+                    wp_cb(ds, (waypoint_t){.depth = depth, .time = waypoint_time, .gas = gas}, SEG_TRAVEL);
+
+                last_waypoint_depth = depth;
+
+                /* switch gas */
+                gas = next;
+
+                ret.tts += add_segment_const(ds, switch_depth, 1, gas);
+
+                if (wp_cb)
+                    wp_cb(ds, (waypoint_t){.depth = depth, .time = 1, .gas = gas}, SEG_GAS_SWITCH);
+
+                continue;
+            }
+
+            /* ascend to current ceiling */
+            ret.tts += add_segment_ascdec(ds, depth, next_stop, fabs(depth - next_stop) / asc_per_min, gas);
+            depth = next_stop;
+            current_gf = get_gf(ds, depth);
+
+            /* if the ceiling moved while we ascended, keep ascending */
+            if (depth > ceiling(ds, current_gf))
+                continue;
+
+            /*
+             * since we've actually reached the ceiling, this is the last of
+             * any number of consecutive travel segments and the waypoint
+             * callback should be called.
+             */
+            waypoint_time = fabs(last_waypoint_depth - depth) / asc_per_min;
+            enum segtype_t segtype = depth <= abs_depth(0) ? SEG_SURFACE : SEG_TRAVEL;
+
+            if (wp_cb)
+                wp_cb(ds, (waypoint_t){.depth = depth, .time = waypoint_time, .gas = gas}, segtype);
+
+            break;
+        }
+
+        /* terminate if we surfaced */
+        if (depth <= abs_depth(0))
+            break;
+
+        /* switch to better gas if available */
+        const gas_t *best = best_gas(depth, deco_gasses, nof_gasses);
+
+        if (best)
+            gas = best;
+
+        /* stop */
+        double stoplen = deco_stop(ds, depth, current_gf, gas);
+        ret.tts += stoplen;
+
+        if (wp_cb)
+            wp_cb(ds, (waypoint_t){.depth = depth, .time = stoplen, .gas = gas}, SEG_DECO_STOP);
+    }
+
+    return ret;
+}
diff --git a/src/schedule.h b/src/schedule.h
new file mode 100644
index 0000000..f2e6363
--- /dev/null
+++ b/src/schedule.h
@@ -0,0 +1,42 @@
+/* SPDX-License-Identifier: MIT-0 */
+
+#ifndef SCHEDULE_H
+#define SCHEDULE_H
+
+#include "deco.h"
+
+typedef struct waypoint_t {
+    double depth;
+    double time;
+    const gas_t *gas;
+} waypoint_t;
+
+typedef struct decoinfo_t {
+    double ndl;
+    double tts;
+} decoinfo_t;
+
+typedef enum segtype_t {
+    SEG_DECO_STOP,
+    SEG_DIVE,
+    SEG_GAS_SWITCH,
+    SEG_NDL,
+    SEG_SAFETY_STOP,
+    SEG_SURFACE,
+    SEG_TRAVEL,
+} segtype_t;
+
+typedef void (*waypoint_callback_t)(const decostate_t *ds, const waypoint_t, const segtype_t);
+
+const gas_t *best_gas(const double depth, const gas_t *gasses, const int nof_gasses);
+const gas_t *next_gas(const double depth, const gas_t *gasses, const int nof_gasses);
+
+int direct_ascent(const decostate_t *ds, const double depth, const double time, const gas_t *gas);
+double calc_ndl(decostate_t *ds, const double depth, const double ascrate, const gas_t *gas);
+
+void simulate_dive(decostate_t *ds, waypoint_t *waypoints, const int nof_waypoints, waypoint_callback_t wp_cb);
+
+decoinfo_t calc_deco(decostate_t *ds, const double start_depth, const gas_t *start_gas, const gas_t *deco_gasses,
+                     const int nof_gasses, waypoint_callback_t wp_cb);
+
+#endif /* end of include guard: SCHEDULE_H */