src/compandt.c

   1 /* libSoX Compander Transfer Function: (c) 2007 robs@users.sourceforge.net
   2  *
   3  * This library is free software; you can redistribute it and/or modify it
   4  * under the terms of the GNU Lesser General Public License as published by
   5  * the Free Software Foundation; either version 2.1 of the License, or (at
   6  * your option) any later version.
   7  *
   8  * This library is distributed in the hope that it will be useful, but
   9  * WITHOUT ANY WARRANTY; without even the implied warranty of
  10  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser
  11  * General Public License for more details.
  12  *
  13  * You should have received a copy of the GNU Lesser General Public License
  14  * along with this library; if not, write to the Free Software Foundation,
  15  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
  16  */
  17
  18 #include "sox_i.h"
  19
  20 #include "compandt.h"
  21 #include <string.h>
  22
  23 #define LOG_TO_LOG10(x) ((x) * 20 / M_LN10)
  24
  25 sox_bool lsx_compandt_show(sox_compandt_t * t, sox_plot_t plot)
  26 {
  27   int i;
  28
  29   for (i = 1; t->segments[i-1].x; ++i)
  30     lsx_debug("TF: %g %g %g %g",
  31        LOG_TO_LOG10(t->segments[i].x),
  32        LOG_TO_LOG10(t->segments[i].y),
  33        LOG_TO_LOG10(t->segments[i].a),
  34        LOG_TO_LOG10(t->segments[i].b));
  35
  36   if (plot == sox_plot_octave) {
  37     printf(
  38       "%% GNU Octave file (may also work with MATLAB(R) )\n"
  39       "in=linspace(-99.5,0,200);\n"
  40       "out=[");
  41     for (i = -199; i <= 0; ++i) {
  42       double in = i/2.;
  43       double in_lin = pow(10., in/20);
  44       printf("%g ", in + 20 * log10(lsx_compandt(t, in_lin)));
  45     }
  46     printf(
  47       "];\n"
  48       "plot(in,out)\n"
  49       "title('SoX effect: compand')\n"
  50       "xlabel('Input level (dB)')\n"
  51       "ylabel('Output level (dB)')\n"
  52       "grid on\n"
  53       "disp('Hit return to continue')\n"
  54       "pause\n");
  55     return sox_false;
  56   }
  57   if (plot == sox_plot_gnuplot) {
  58     printf(
  59       "# gnuplot file\n"
  60       "set title 'SoX effect: compand'\n"
  61       "set xlabel 'Input level (dB)'\n"
  62       "set ylabel 'Output level (dB)'\n"
  63       "set grid xtics ytics\n"
  64       "set key off\n"
  65       "plot '-' with lines\n");
  66     for (i = -199; i <= 0; ++i) {
  67       double in = i/2.;
  68       double in_lin = pow(10., in/20);
  69       printf("%g %g\n", in, in + 20 * log10(lsx_compandt(t, in_lin)));
  70     }
  71     printf(
  72       "e\n"
  73       "pause -1 'Hit return to continue'\n");
  74     return sox_false;
  75   }
  76   return sox_true;
  77 }
  78
  79 static void prepare_transfer_fn(sox_compandt_t * t)
  80 {
  81   int i;
  82   double radius = t->curve_dB * M_LN10 / 20;
  83
  84   for (i = 0; !i || t->segments[i-2].x; i += 2) {
  85     t->segments[i].y += t->outgain_dB;
  86     t->segments[i].x *= M_LN10 / 20; /* Convert to natural logs */
  87     t->segments[i].y *= M_LN10 / 20;
  88   }
  89
  90 #define line1 t->segments[i - 4]
  91 #define curve t->segments[i - 3]
  92 #define line2 t->segments[i - 2]
  93 #define line3 t->segments[i - 0]
  94   for (i = 4; t->segments[i - 2].x; i += 2) {
  95     double x, y, cx, cy, in1, in2, out1, out2, theta, len, r;
  96
  97     line1.a = 0;
  98     line1.b = (line2.y - line1.y) / (line2.x - line1.x);
  99
 100     line2.a = 0;
 101     line2.b = (line3.y - line2.y) / (line3.x - line2.x);
 102
 103     theta = atan2(line2.y - line1.y, line2.x - line1.x);
 104     len = sqrt(pow(line2.x - line1.x, 2.) + pow(line2.y - line1.y, 2.));
 105     r = min(radius, len);
 106     curve.x = line2.x - r * cos(theta);
 107     curve.y = line2.y - r * sin(theta);
 108
 109     theta = atan2(line3.y - line2.y, line3.x - line2.x);
 110     len = sqrt(pow(line3.x - line2.x, 2.) + pow(line3.y - line2.y, 2.));
 111     r = min(radius, len / 2);
 112     x = line2.x + r * cos(theta);
 113     y = line2.y + r * sin(theta);
 114
 115     cx = (curve.x + line2.x + x) / 3;
 116     cy = (curve.y + line2.y + y) / 3;
 117
 118     line2.x = x;
 119     line2.y = y;
 120
 121     in1 = cx - curve.x;
 122     out1 = cy - curve.y;
 123     in2 = line2.x - curve.x;
 124     out2 = line2.y - curve.y;
 125     curve.a = (out2/in2 - out1/in1) / (in2-in1);
 126     curve.b = out1/in1 - curve.a*in1;
 127   }
 128 #undef line1
 129 #undef curve
 130 #undef line2
 131 #undef line3
 132   t->segments[i - 3].x = 0;
 133   t->segments[i - 3].y = t->segments[i - 2].y;
 134
 135   t->in_min_lin = exp(t->segments[1].x);
 136   t->out_min_lin= exp(t->segments[1].y);
 137 }
 138
 139 static sox_bool parse_transfer_value(char const * text, double * value)
 140 {
 141   char dummy;     /* To check for extraneous chars. */
 142
 143   if (!text) {
 144     lsx_fail("syntax error trying to read transfer function value");
 145     return sox_false;
 146   }
 147   if (!strcmp(text, "-inf"))
 148     *value = -20 * log10(-(double)SOX_SAMPLE_MIN);
 149   else if (sscanf(text, "%lf %c", value, &dummy) != 1) {
 150     lsx_fail("syntax error trying to read transfer function value");
 151     return sox_false;
 152   }
 153   else if (*value > 0) {
 154     lsx_fail("transfer function values are relative to maximum volume so can't exceed 0dB");
 155     return sox_false;
 156   }
 157   return sox_true;
 158 }
 159
 160 sox_bool lsx_compandt_parse(sox_compandt_t * t, char * points, char * gain)
 161 {
 162   char const * text = points;
 163   unsigned i, j, num, pairs, commas = 0;
 164   char dummy;     /* To check for extraneous chars. */
 165
 166   if (sscanf(points, "%lf %c", &t->curve_dB, &dummy) == 2 && dummy == ':')
 167     points = strchr(points, ':') + 1;
 168   else t->curve_dB = 0;
 169   t->curve_dB = max(t->curve_dB, .01);
 170
 171   while (*text) commas += *text++ == ',';
 172   pairs = 1 + commas / 2;
 173   ++pairs;    /* allow room for extra pair at the beginning */
 174   pairs *= 2; /* allow room for the auto-curves */
 175   ++pairs;    /* allow room for 0,0 at end */
 176   t->segments = lsx_calloc(pairs, sizeof(*t->segments));
 177
 178 #define s(n) t->segments[2*((n)+1)]
 179   for (i = 0, text = strtok(points, ","); text != NULL; ++i) {
 180     if (!parse_transfer_value(text, &s(i).x))
 181       return sox_false;
 182     if (i && s(i-1).x > s(i).x) {
 183       lsx_fail("transfer function input values must be strictly increasing");
 184       return sox_false;
 185     }
 186     if (i || (commas & 1)) {
 187       text = strtok(NULL, ",");
 188       if (!parse_transfer_value(text, &s(i).y))
 189         return sox_false;
 190       s(i).y -= s(i).x;
 191     }
 192     text = strtok(NULL, ",");
 193   }
 194   num = i;
 195
 196   if (num == 0 || s(num-1).x) /* Add 0,0 if necessary */
 197     ++num;
 198 #undef s
 199
 200   if (gain && sscanf(gain, "%lf %c", &t->outgain_dB, &dummy) != 1) {
 201     lsx_fail("syntax error trying to read post-processing gain value");
 202     return sox_false;
 203   }
 204
 205 #define s(n) t->segments[2*(n)]
 206   s(0).x = s(1).x - 2 * t->curve_dB; /* Add a tail off segment at the start */
 207   s(0).y = s(1).y;
 208   ++num;
 209
 210   for (i = 2; i < num; ++i) { /* Join adjacent colinear segments */
 211     double g1 = (s(i-1).y - s(i-2).y) * (s(i-0).x - s(i-1).x);
 212     double g2 = (s(i-0).y - s(i-1).y) * (s(i-1).x - s(i-2).x);
 213     if (fabs(g1 - g2)) /* fabs stops epsilon problems */
 214       continue;
 215     --num;
 216     for (j = --i; j < num; ++j)
 217       s(j) = s(j+1);
 218   }
 219 #undef s
 220
 221   prepare_transfer_fn(t);
 222   return sox_true;
 223 }
 224
 225 void lsx_compandt_kill(sox_compandt_t * p)
 226 {
 227   free(p->segments);
 228 }
 229