diff options
author | Dirk Hohndel <dirk@hohndel.org> | 2012-11-05 08:56:18 -0800 |
---|---|---|
committer | Dirk Hohndel <dirk@hohndel.org> | 2012-11-05 08:56:18 -0800 |
commit | 853277ba9d821922a8da6e55dd47d2f09def50fa (patch) | |
tree | 31b22bc4c815acc910319595bc61ba46bbaa05b2 /profile.c | |
parent | 5f2f415cdbc7854fd6554a75bcf37d9462747a00 (diff) | |
download | subsurface-853277ba9d821922a8da6e55dd47d2f09def50fa.tar.gz |
Plot text values for partial pressure graphs
The algorithms attempt to identify "interesting" points where the user
might want to know the value of the graph.
Signed-off-by: Dirk Hohndel <dirk@hohndel.org>
Diffstat (limited to 'profile.c')
-rw-r--r-- | profile.c | 263 |
1 files changed, 262 insertions, 1 deletions
@@ -8,6 +8,7 @@ #include <stdarg.h> #include <string.h> #include <time.h> +#include <math.h> #include "dive.h" #include "display.h" @@ -505,6 +506,264 @@ static void plot_depth_scale(struct graphics_context *gc, struct plot_info *pi) } } +/* ap points to an array of int with pi->nr + 1 elements that is + * ininitialized with just one -1 entry + * this adds entries (if they aren't too close to an existing one) + * and keeps things sorted + * we KNOW the array is big enough to hold all possible indices + * a2p is a secondary array - we insert value at the same relative + * positio as idx in ap */ +static void add_index(int idx, int margin, int **ap, int **a2p, int value) +{ + int j, i = 0; + int *a = *ap; + int *a2 = *a2p; + + while (a[i] != -1 && a[i] < idx) + i++; + if (a[i] == idx) + return; + if (a[i] != -1 && a[i - 1] != -1 && idx - a[i - 1] < margin) + return; + if (a[i] != -1 && a[i] - idx < margin) + return; + j = i; + while (a[j] != -1) + j++; + while (j >= i) { + a[j+1] = a[j]; + a2[j+1] = a2[j]; + j--; + } + a[i] = idx; + a2[i] = value; +} + +#define LI(_i,_j) MAX((_i)-(_j), 0) +#define RI(_i,_j) MIN((_i)+(_j), nr - 1) +#define SPIKE(_i,_s) if (fabs(_s) > fabs(spk_data[_i])) spk_data[_i] = (_s) +/* this is an attempt at a metric that finds spikes in a data series */ +static void calculate_spikyness(int nr, double *data, double *spk_data, int deltax, double deltay) +{ + int i, j; + double dminl, dminr, dmaxl, dmaxr; + +#if DEBUG_PROFILE > 2 + printf("Spike data: \n 0 "); +#endif + for (i = 0; i < nr; i++) { + dminl = dminr = dmaxl = dmaxr = data[i]; + spk_data[i] = 0.0; + for (j = 1; j < deltax; j++) { + if (data[LI(i,j)] < dminl) + dminl = data[LI(i,j)]; + if (data[LI(i,j)] > dmaxl) + dmaxl = data[LI(i,j)]; + + if (data[RI(i,j)] < dminr) + dminr = data[RI(i,j)]; + if (data[RI(i,j)] > dmaxr) + dmaxr = data[RI(i,j)]; + + /* don't do super narrow */ + if (j < deltax / 3) + continue; + /* falling edge on left */ + if (dmaxl == data[i] && dmaxr - data[i] < 0.1 * (data[i] - dminl)) + SPIKE(i,(data[i] - dminl) / j); + /* falling edge on right */ + if (dmaxr == data[i] && dmaxl - data[i] < 0.1 * (data[i] - dminr)) + SPIKE(i,(data[i] - dminr) / j); + + /* minima get a negative spike value */ + /* rising edge on left */ + if (dminl == data[i] && data[i] - dminr < 0.1 * (dmaxl - data[i])) + SPIKE(i,(data[i] - dmaxl) / j); + /* rising edge on right */ + if (dminr == data[i] && data[i] - dminl < 0.1 * (dmaxr - data[i])) + SPIKE(i,(data[i] - dmaxr) / j); + } +#if DEBUG_PROFILE > 2 + fprintf(debugfile, "%.4lf ", spk_data[i]); + if (i % 12 == 11) + fprintf(debugfile, "\n%2d ", (i + 1) / 12); +#endif + } +#if DEBUG_PROFILE > 2 + printf("\n"); +#endif +} + +/* only show one spike in a deltax wide region - pick the highest (and first if the same) */ +static gboolean higher_spike(double *spk_data, int idx, int nr, int deltax) +{ + int i; + double s = fabs(spk_data[idx]); + for (i = MAX(0, idx - deltax); i <= MIN(idx + deltax, nr - 1); i++) + if (fabs(spk_data[i]) > s) + return TRUE; + else if (fabs(spk_data[i]) == s && i < idx) + return TRUE; + return FALSE; +} + +/* this figures out which time stamps provide "interesting" formations in the graphs; + * this includes local minima and maxima as well as long plateaus. + * pass in the function that returns the value at a certain point (as double), + * the delta in time (expressed as number of data points of "significant time") + * the delta at which the value is considered to have been "significantly changed" and + * the number of points to cover + * returns a list of indices that ends with a -1 of times that are "interesting" */ +static void find_points_of_interest(struct plot_info *pi, double (*value_func)(int, struct plot_info *), + int deltax, double deltay, int **poip, int **poip_vpos) +{ + int i, j, nr = pi->nr; + double *data, *data_max, *data_min, *spk_data; + double min, max; + int *pois; + + /* avoid all the function calls by creating a local array and + * have some helper arrays to make our lifes easier */ + + data = malloc(nr * sizeof(double)); + data_max = malloc(nr * sizeof(double)); + data_min = malloc(nr * sizeof(double)); + spk_data = malloc(nr * sizeof(double)); + + pois = *poip = malloc((nr + 1) * sizeof(int)); + *poip_vpos = malloc((nr + 1) * sizeof(int)); + pois[0] = -1; + pois[1] = -1; + + /* copy the data and get the absolute minimum and maximum while we do it */ + for (i = 0; i < nr; i++) { + data_max[i] = data_min[i] = data[i] = value_func(i, pi); + if (i == 0 || data[i] < min) + min = data[i]; + if (i == 0 || data[i] > max) + max = data[i]; + } + /* next find out if there are real spikes in the graph */ + calculate_spikyness(nr, data, spk_data, deltax, deltay); + + /* now process all data points */ + for (i = 0; i < nr; i++) { + /* get the local min/max */ + for (j = MAX(0, i - deltax); j < i + deltax && j < nr; j++) { + if (data[j] < data[i]) + data_min[i] = data[j]; + if (data[j] > data[i]) + data_max[i] = data[j]; + } + /* is i the overall minimum or maximum */ + if (data[i] == max) + add_index(i, deltax, poip, poip_vpos, BOTTOM); + if (data[i] == min) + add_index(i, deltax, poip, poip_vpos, TOP); + /* is i a spike? */ + if (fabs(spk_data[i]) > 0.01 && ! higher_spike(spk_data, i, nr, deltax)) { + if (spk_data[i] > 0.0) + add_index(i, deltax, poip, poip_vpos, BOTTOM); + if (spk_data[i] < 0.0) + add_index(i, deltax, poip, poip_vpos, TOP); + } + /* is i a significant local minimum or maximum? */ + if (data[i] == data_min[i] && data_max[i] - data[i] > deltay) + add_index(i, deltax, poip, poip_vpos, TOP); + if (data[i] == data_max[i] && data[i] - data_min[i] > deltay) + add_index(i, deltax, poip, poip_vpos, BOTTOM); + } + /* still need to search for plateaus */ +} + +static void setup_pp_limits(struct graphics_context *gc, struct plot_info *pi) +{ + int maxdepth; + + gc->leftx = 0; + gc->rightx = get_maxtime(pi); + + /* the maxdepth already includes extra vertical space - and if + * we take the corresponding pressure as maximum partial + * pressure the graph seems to look fine*/ + maxdepth = get_maxdepth(pi); + gc->topy = (maxdepth + 10000) / 10000.0 * 1.01325; + gc->bottomy = 0.0; +} + +static void plot_single_pp_text(struct graphics_context *gc, int sec, double pp, color_indice_t color) +{ + text_render_options_t tro = {12, color, CENTER, BOTTOM}; + plot_text(gc, &tro, sec, pp, "%.1lf", pp); + + if (color == PN2) + printf("pN2 %lf\n", pp); +} + +#define MAXPP(_mpp, _pp) { _mpp = 0; \ + for(i = 0; i< pi->nr; i++) \ + if (pi->entry[i]._pp > _mpp) \ + _mpp = pi->entry[i]._pp; \ + } + +static double po2_value(int idx, struct plot_info *pi) +{ + return pi->entry[idx].po2; +} + +static double pn2_value(int idx, struct plot_info *pi) +{ + return pi->entry[idx].pn2; +} + +static double phe_value(int idx, struct plot_info *pi) +{ + return pi->entry[idx].phe; +} + +static void plot_single_gas_pp_text(struct graphics_context *gc, struct plot_info *pi, + double (*value_func)(int, struct plot_info *), + double value_threshold, int color) +{ + int *pois, *pois_vpos; + int i, two_minutes = 1; + /* don't bother with local min/max if the dive is under two minutes */ + if (pi->entry[pi->nr - 1].sec > 120) { + int idx = 0; + while (pi->entry[idx].sec == 0) + idx++; + while (pi->entry[idx + two_minutes].sec < 120) + two_minutes++; + } else { + two_minutes = pi->nr; + } + find_points_of_interest(pi, value_func, two_minutes, value_threshold, &pois, &pois_vpos); + for (i = 0; pois[i] != -1; i++) { + struct plot_data *entry = pi->entry + pois[i]; +#if DEBUG_PROFILE > 1 + fprintf(debugfile, "POI at %d sec value %lf\n", entry->sec, entry->po2); +#endif + plot_single_pp_text(gc, entry->sec, value_func(pois[i], pi), color); + } + free(pois); + free(pois_vpos); +} + +static void plot_pp_text(struct graphics_context *gc, struct plot_info *pi) +{ + setup_pp_limits(gc, pi); + + if (enabled_graphs.po2) { + plot_single_gas_pp_text(gc, pi, po2_value, 0.4, PO2); + } + if (enabled_graphs.pn2) { + plot_single_gas_pp_text(gc, pi, pn2_value, 0.4, PN2); + } + if (enabled_graphs.phe) { + plot_single_gas_pp_text(gc, pi, phe_value, 0.4, PHE); + } +} + static void plot_pp_gas_profile(struct graphics_context *gc, struct plot_info *pi) { int i; @@ -1591,8 +1850,10 @@ void plot(struct graphics_context *gc, cairo_rectangle_t *drawing_area, struct d cairo_close_path(gc->cr); cairo_stroke(gc->cr); - if (GRAPHS_ENABLED) + if (GRAPHS_ENABLED) { plot_pp_gas_profile(gc, pi); + plot_pp_text(gc, pi); + } /* now shift the translation back by half the margin; * this way we can draw the vertical scales on both sides */ |