diff options
-rw-r--r-- | deco.c | 63 | ||||
-rw-r--r-- | dive.h | 5 | ||||
-rw-r--r-- | divelist.c | 77 | ||||
-rw-r--r-- | divelist.h | 1 | ||||
-rw-r--r-- | profile.c | 9 |
5 files changed, 137 insertions, 18 deletions
@@ -7,10 +7,11 @@ * The implementation below is (C) Dirk Hohndel 2012 and released under the GPLv2 * * clear_deco() - call to initialize for a new deco calculation - * add_segment(pressure, gasmix) - add 1 second at the given pressure, breathing gasmix + * add_segment(pressure, gasmix, seconds) - add <seconds> at the given pressure, breathing gasmix * deco_allowed_depth(tissues_tolerance, surface_pressure, dive, smooth) * - ceiling based on lead tissue, surface pressure, 3m increments or smooth */ +#include <math.h> #include "dive.h" //! Option structure for Buehlmann decompression. @@ -76,6 +77,7 @@ const double buehlmann_He_factor_expositon_one_second[] = { 1.00198406028040E-004, 7.83611475491108E-005, 6.13689891868496E-005, 4.81280465299827E-005}; #define WV_PRESSURE 0.0627 /* water vapor pressure */ +#define N2_IN_AIR 0.7902 #define DIST_FROM_3_MTR 0.28 #define PRESSURE_CHANGE_3M 0.3 #define TOLERANCE 0.02 @@ -145,32 +147,65 @@ static double tissue_tolerance_calc(void) } /* add a second at the given pressure and gas to the deco calculation */ -double add_segment(double pressure, struct gasmix *gasmix) +double add_segment(double pressure, struct gasmix *gasmix, int period_in_seconds) { int ci; double ppn2 = (pressure - WV_PRESSURE) * (1000 - gasmix->o2.permille - gasmix->he.permille) / 1000.0; double pphe = (pressure - WV_PRESSURE) * gasmix->he.permille / 1000.0; - divetime++; /* right now we just do OC */ - for (ci = 0; ci < 16; ci++) { - if (ppn2 - tissue_n2_sat[ci] > 0) - tissue_n2_sat[ci] += buehlmann_config.satmult * (ppn2 - tissue_n2_sat[ci]) * buehlmann_N2_factor_expositon_one_second[ci]; - else - tissue_n2_sat[ci] += buehlmann_config.desatmult * (ppn2 - tissue_n2_sat[ci]) * buehlmann_N2_factor_expositon_one_second[ci]; - if (pphe - tissue_he_sat[ci] > 0) - tissue_he_sat[ci] += buehlmann_config.satmult * (pphe - tissue_he_sat[ci]) * buehlmann_He_factor_expositon_one_second[ci]; - else - tissue_he_sat[ci] += buehlmann_config.desatmult * (pphe - tissue_he_sat[ci]) * buehlmann_He_factor_expositon_one_second[ci]; + if (period_in_seconds == 1) { /* that's what we do during the dive */ + for (ci = 0; ci < 16; ci++) { + if (ppn2 - tissue_n2_sat[ci] > 0) + tissue_n2_sat[ci] += buehlmann_config.satmult * (ppn2 - tissue_n2_sat[ci]) * + buehlmann_N2_factor_expositon_one_second[ci]; + else + tissue_n2_sat[ci] += buehlmann_config.desatmult * (ppn2 - tissue_n2_sat[ci]) * + buehlmann_N2_factor_expositon_one_second[ci]; + if (pphe - tissue_he_sat[ci] > 0) + tissue_he_sat[ci] += buehlmann_config.satmult * (pphe - tissue_he_sat[ci]) * + buehlmann_He_factor_expositon_one_second[ci]; + else + tissue_he_sat[ci] += buehlmann_config.desatmult * (pphe - tissue_he_sat[ci]) * + buehlmann_He_factor_expositon_one_second[ci]; + } + } else { /* all other durations */ + for (ci = 0; ci < 16; ci++) + { + if (ppn2 - tissue_n2_sat[ci] > 0) + tissue_n2_sat[ci] += buehlmann_config.satmult * (ppn2 - tissue_n2_sat[ci]) * + (1 - pow(2.0,(- period_in_seconds / (buehlmann_N2_t_halflife[ci] * 60)))); + else + tissue_n2_sat[ci] += buehlmann_config.desatmult * (ppn2 - tissue_n2_sat[ci]) * + (1 - pow(2.0,(- period_in_seconds / (buehlmann_N2_t_halflife[ci] * 60)))); + if (pphe - tissue_he_sat[ci] > 0) + tissue_he_sat[ci] += buehlmann_config.satmult * (pphe - tissue_he_sat[ci]) * + (1 - pow(2.0,(- period_in_seconds / (buehlmann_He_t_halflife[ci] * 60)))); + else + tissue_he_sat[ci] += buehlmann_config.desatmult * (pphe - tissue_he_sat[ci]) * + (1 - pow(2.0,(- period_in_seconds / (buehlmann_He_t_halflife[ci] * 60)))); + } } return tissue_tolerance_calc(); } -void clear_deco() +void dump_tissues() +{ + int ci; + printf("N2 tissues:"); + for (ci = 0; ci < 16; ci++) + printf(" %6.3e", tissue_n2_sat[ci]); + printf("\nHe tissues:"); + for (ci = 0; ci < 16; ci++) + printf(" %6.3e", tissue_he_sat[ci]); + printf("\n"); +} + +void clear_deco(double surface_pressure) { int ci; for (ci = 0; ci < 16; ci++) { - tissue_n2_sat[ci] = 0.0; + tissue_n2_sat[ci] = (surface_pressure - WV_PRESSURE) * N2_IN_AIR; tissue_he_sat[ci] = 0.0; tissue_tolerated_ambient_pressure[ci] = 0.0; } @@ -572,8 +572,9 @@ extern void subsurface_command_line_exit(gint *, gchar ***); #define FRACTION(n,x) ((unsigned)(n)/(x)),((unsigned)(n)%(x)) -extern double add_segment(double pressure, struct gasmix *gasmix); -extern void clear_deco(void); +extern double add_segment(double pressure, struct gasmix *gasmix, int period_in_seconds); +extern void clear_deco(double surface_pressure); +extern void dump_tissues(void); extern unsigned int deco_allowed_depth(double tissues_tolerance, double surface_pressure, struct dive *dive, gboolean smooth); #ifdef DEBUGFILE extern char *debugfilename; diff --git a/divelist.c b/divelist.c index 5df8fbb57..32f5f44a3 100644 --- a/divelist.c +++ b/divelist.c @@ -818,6 +818,83 @@ static int calculate_sac(struct dive *dive, struct divecomputer *dc) return sac * 1000; } +/* for now we do this based on the first divecomputer */ +static void add_dive_to_deco(struct dive *dive) +{ + struct divecomputer *dc = &dive->dc; + int i; + + if (!dc) + return; + for (i = 1; i < dive->dc.samples; i++) { + struct sample *psample = dc->sample + i - 1; + struct sample *sample = dc->sample + i; + int t0 = psample->time.seconds; + int t1 = sample->time.seconds; + int j; + + for (j = t0; j < t1; j++) { + int depth = 0.5 + psample->depth.mm + (j - t0) * (sample->depth.mm - psample->depth.mm) / (t1 - t0); + (void) add_segment(depth_to_mbar(depth, dive) / 1000.0, &dive->cylinder[sample->sensor].gasmix, 1); + } + } +} + +static struct gasmix air = { .o2.permille = 209 }; + +/* take into account previous dives until there is a 48h gap between dives */ +void init_decompression(struct dive *dive) +{ + int i, divenr = -1; + timestamp_t when; + gboolean deco_init = FALSE; + + if (!dive) + return; + while (++divenr < dive_table.nr && get_dive(divenr) != dive) + ; + when = dive->when; + i = divenr; + while (--i) { + struct dive* pdive = get_dive(i); + if (!pdive || pdive->when > when || pdive->when + pdive->duration.seconds + 48 * 60 * 60 < when) + break; + when = pdive->when; + } + + while (++i < divenr) { + struct dive* pdive = get_dive(i); + double surface_pressure = pdive->surface_pressure.mbar ? pdive->surface_pressure.mbar / 1000.0 : 1.013; + unsigned int surface_time = get_dive(i+1)->when - pdive->when - pdive->duration.seconds; + + if (!deco_init) { + clear_deco(surface_pressure); + deco_init = TRUE; +#if DEBUG & 16 + dump_tissues(); +#endif + } + add_dive_to_deco(pdive); +#if DEBUG & 16 + printf("added dive #%d\n", pdive->number); + dump_tissues(); +#endif + add_segment(surface_pressure, &air, surface_time); +#if DEBUG & 16 + printf("after surface intervall of %d:%02u\n", FRACTION(surface_time,60)); + dump_tissues(); +#endif + } + if (!deco_init) { + double surface_pressure = dive->surface_pressure.mbar ? dive->surface_pressure.mbar / 1000.0 : 1.013; + clear_deco(surface_pressure); +#if DEBUG & 16 + printf("no previous dive\n"); + dump_tissues(); +#endif + } +} + void update_cylinder_related_info(struct dive *dive) { if (dive != NULL) { diff --git a/divelist.h b/divelist.h index b309e9a54..41a9f38ae 100644 --- a/divelist.h +++ b/divelist.h @@ -15,4 +15,5 @@ extern void remember_tree_state(void); extern void restore_tree_state(void); extern void select_next_dive(void); extern void select_prev_dive(void); +extern void init_decompression(struct dive * dive); #endif @@ -1566,7 +1566,8 @@ static struct plot_info *create_plot_info(struct dive *dive, struct divecomputer pi = &gc->pi; /* reset deco information to start the calculation */ - clear_deco(); + init_decompression(dive); + /* we want to potentially add synthetic plot_info elements for the gas changes */ nr = dc->samples + 4 + 2 * count_gas_change_events(dc); if (last_pi_entry) @@ -1756,13 +1757,17 @@ static struct plot_info *create_plot_info(struct dive *dive, struct divecomputer float ceiling_pressure = 0; for (j = t0; j < t1; j++) { int depth = 0.5 + (entry - 1)->depth + (j - t0) * (entry->depth - (entry - 1)->depth) / (t1 - t0); - double min_pressure = add_segment(depth_to_mbar(depth, dive) / 1000.0, &dive->cylinder[cylinderindex].gasmix); + double min_pressure = add_segment(depth_to_mbar(depth, dive) / 1000.0, + &dive->cylinder[cylinderindex].gasmix, 1); if (min_pressure > ceiling_pressure) ceiling_pressure = min_pressure; } entry->ceiling = deco_allowed_depth(ceiling_pressure, surface_pressure, dive, !prefs.calc_ceiling_3m_incr); } } +#if DECO_CALC_DEBUG + dump_tissues(); +#endif if (entry) current->t_end = entry->sec; |