cinelerra-5.0/cinelerra/resample.C

   1
   2 /*
   3  * CINELERRA
   4  * Copyright (C) 2009 Adam Williams <broadcast at earthling dot net>
   5  *
   6  * This program is free software; you can redistribute it and/or modify
   7  * it under the terms of the GNU General Public License as published by
   8  * the Free Software Foundation; either version 2 of the License, or
   9  * (at your option) any later version.
  10  *
  11  * This program is distributed in the hope that it will be useful,
  12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14  * GNU General Public License for more details.
  15  *
  16  * You should have received a copy of the GNU General Public License
  17  * along with this program; if not, write to the Free Software
  18  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  19  *
  20  */
  21
  22 #include "bcsignals.h"
  23 #include "clip.h"
  24 #include "format.inc"
  25 #include "resample.h"
  26 #include "samples.h"
  27 #include "transportque.inc"
  28
  29 #include <math.h>
  30 #include <stdio.h>
  31 #include <stdlib.h>
  32 #include <string.h>
  33
  34 // Resampling from Lame
  35
  36 Resample::Resample()
  37 {
  38         old = new double[BLACKSIZE];
  39         resample_init = 0;
  40         last_ratio = 0;
  41         output_temp = 0;
  42         output_size = 0;
  43         output_allocation = 0;
  44         input_size = RESAMPLE_CHUNKSIZE;
  45         input_position = 0;
  46         input = new Samples(input_size + 1);
  47         output_position = 0;
  48         itime = 0;
  49         direction = PLAY_FORWARD;
  50 }
  51
  52
  53 Resample::~Resample()
  54 {
  55         delete [] old;
  56         delete [] output_temp;
  57         delete input;
  58 }
  59
  60 int Resample::read_samples(Samples *buffer, int64_t start, int64_t len)
  61 {
  62         return 0;
  63 }
  64
  65 int Resample::get_direction()
  66 {
  67         return direction;
  68 }
  69
  70
  71 void Resample::reset()
  72 {
  73         resample_init = 0;
  74         output_size = 0;
  75         output_position = 0;
  76         input_position = 0;
  77 }
  78
  79 double Resample::blackman(int i, double offset, double fcn, int l)
  80 {
  81   /* This algorithm from:
  82 SIGNAL PROCESSING ALGORITHMS IN FORTRAN AND C
  83 S.D. Stearns and R.A. David, Prentice-Hall, 1992
  84   */
  85
  86         double bkwn;
  87         double wcn = (M_PI * fcn);
  88         double dly = l / 2.0;
  89         double x = i-offset;
  90         if(x < 0) x = 0;
  91         else
  92         if(x > l) x = l;
  93
  94         bkwn = 0.42 - 0.5 * cos((x * 2) * M_PI /l) + 0.08 * cos((x * 4) * M_PI /l);
  95         if(fabs(x - dly) < 1e-9)
  96                 return wcn / M_PI;
  97     else
  98         return (sin((wcn * (x - dly))) / (M_PI * (x - dly)) * bkwn);
  99 }
 100
 101
 102 int Resample::get_output_size()
 103 {
 104         return output_size;
 105 }
 106
 107 // void Resample::read_output(double *output, int size)
 108 // {
 109 //      memcpy(output, output_temp, size * sizeof(double));
 110 // // Shift leftover forward
 111 //      for(int i = size; i < output_size; i++)
 112 //              output_temp[i - size] = output_temp[i];
 113 //      output_size -= size;
 114 // }
 115
 116
 117
 118 void Resample::resample_chunk(Samples *input_buffer,
 119         int64_t in_len,
 120         int in_rate,
 121         int out_rate)
 122 {
 123         double resample_ratio = (double)in_rate / out_rate;
 124         int filter_l;
 125         double fcn, intratio;
 126         double offset, xvalue;
 127         int num_used;
 128         int i, j, k;
 129         double *input = input_buffer->get_data();
 130 //printf("Resample::resample_chunk %d in_len=" _LD " input_size=%d\n",
 131 //__LINE__,
 132 //in_len,
 133 //input_size);
 134
 135         intratio = (fabs(resample_ratio - floor(.5 + resample_ratio)) < .0001);
 136         fcn = .90 / resample_ratio;
 137         if(fcn > .90) fcn = .90;
 138         filter_l = BLACKSIZE - 6;
 139 /* must be odd */
 140         if(0 == filter_l % 2 ) --filter_l;
 141
 142 /* if resample_ratio = int, filter_l should be even */
 143         filter_l += (int)intratio;
 144
 145 // Blackman filter initialization must be called whenever there is a
 146 // sampling ratio change
 147         if(!resample_init || last_ratio != resample_ratio)
 148         {
 149                 resample_init = 1;
 150                 itime = 0;
 151                 bzero(old, sizeof(double) * BLACKSIZE);
 152
 153 // precompute blackman filter coefficients
 154         for (j = 0; j <= 2 * BPC; ++j)
 155                 {
 156                         for(j = 0; j <= 2 * BPC; j++)
 157                         {
 158                                 offset = (double)(j - BPC) / (2 * BPC);
 159                                 for(i = 0; i <= filter_l; i++)
 160                                 {
 161                                         blackfilt[j][i] = blackman(i, offset, fcn, filter_l);
 162                                 }
 163                         }
 164                 }
 165         }
 166
 167 // Main loop
 168         double *inbuf_old = old;
 169         for(k = 0; 1; k++)
 170         {
 171                 double time0;
 172                 int joff;
 173
 174                 time0 = k * resample_ratio;
 175                 j = (int)floor(time0 - itime);
 176
 177 //              if(j + filter_l / 2 >= input_size) break;
 178                 if(j + filter_l / 2 >= in_len) break;
 179
 180 /* blackman filter.  by default, window centered at j+.5(filter_l%2) */
 181 /* but we want a window centered at time0.   */
 182                 offset = (time0 - itime - (j + .5 * (filter_l % 2)));
 183                 joff = (int)floor((offset * 2 * BPC) + BPC + .5);
 184                 xvalue = 0;
 185
 186
 187                 for(i = 0; i <= filter_l; i++)
 188                 {
 189                         int j2 = i + j - filter_l / 2;
 190 //printf("j2=%d\n", j2);
 191                         double y = ((j2 < 0) ? inbuf_old[BLACKSIZE + j2] : input[j2]);
 192
 193                         xvalue += y * blackfilt[joff][i];
 194                 }
 195
 196
 197                 if(output_allocation <= output_size)
 198                 {
 199                         double *new_output = 0;
 200                         int64_t new_allocation = output_allocation ? (output_allocation * 2) : 16384;
 201                         new_output = new double[new_allocation];
 202                         if(output_temp)
 203                         {
 204                                 bcopy(output_temp, new_output, output_allocation * sizeof(double));
 205                                 delete [] output_temp;
 206                         }
 207
 208                         output_temp = new_output;
 209                         output_allocation = new_allocation;
 210                 }
 211
 212                 output_temp[output_size++] = xvalue;
 213         }
 214
 215         num_used = MIN(in_len, j + filter_l / 2);
 216         itime += num_used - k * resample_ratio;
 217         for(i = 0; i < BLACKSIZE; i++)
 218                 inbuf_old[i] = input[num_used + i - BLACKSIZE];
 219
 220         last_ratio = resample_ratio;
 221
 222 }
 223
 224 int Resample::read_chunk(Samples *input,
 225         int64_t len)
 226 {
 227         int fragment = len;
 228         if(direction == PLAY_REVERSE &&
 229                 input_position - len < 0)
 230         {
 231                 fragment = input_position;
 232         }
 233
 234         int result = read_samples(input, input_position, fragment);
 235
 236         if(direction == PLAY_FORWARD)
 237         {
 238                 input_position += fragment;
 239         }
 240         else
 241         {
 242                 input_position -= fragment;
 243 // Mute unused part of buffer
 244                 if(fragment < len)
 245                 {
 246                         bzero(input->get_data() + fragment,
 247                                 (len - fragment) * sizeof(double));
 248                 }
 249         }
 250
 251         return result;
 252 }
 253
 254
 255 void Resample::reverse_buffer(double *buffer, int64_t len)
 256 {
 257         double *ap = buffer;
 258         double *bp = ap + len;
 259         while( ap < --bp ) {
 260                 double t = *ap;
 261                 *ap++ = *bp;
 262                 *bp = t;
 263         }
 264 }
 265
 266
 267 int Resample::resample(Samples *output,
 268         int64_t out_len,
 269         int in_rate,
 270         int out_rate,
 271         int64_t out_position,
 272         int direction)
 273 {
 274         int result = 0;
 275
 276
 277 // printf("Resample::resample 1 output_position=" _LD " out_position=" _LD " out_len=" _LD "\n",
 278 // output_position,
 279 // out_position,
 280 // out_len);
 281 // Changed position
 282         if(labs(this->output_position - out_position) > 0 ||
 283                 direction != this->direction)
 284         {
 285                 reset();
 286
 287 // Compute starting point in input rate.
 288                 this->input_position = out_position * in_rate / out_rate;
 289                 this->direction = direction;
 290         }
 291
 292
 293         int remaining_len = out_len;
 294         double *output_ptr = output->get_data();
 295         while(remaining_len > 0 && !result)
 296         {
 297 // Drain output buffer
 298                 if(output_size)
 299                 {
 300                         int fragment_len = output_size;
 301                         if(fragment_len > remaining_len) fragment_len = remaining_len;
 302
 303 //printf("Resample::resample 1 %d %d\n", remaining_len, output_size);
 304                         bcopy(output_temp, output_ptr, fragment_len * sizeof(double));
 305
 306 // Shift leftover forward
 307                         for(int i = fragment_len; i < output_size; i++)
 308                                 output_temp[i - fragment_len] = output_temp[i];
 309
 310                         output_size -= fragment_len;
 311                         remaining_len -= fragment_len;
 312                         output_ptr += fragment_len;
 313                 }
 314
 315 // Import new samples
 316 //printf("Resample::resample 2 %d\n", remaining_len);
 317                 if(remaining_len > 0)
 318                 {
 319 //printf("Resample::resample 3 input_size=%d out_position=%d\n", input_size, out_position);
 320                         result = read_chunk(input, input_size);
 321                         resample_chunk(input,
 322                                 input_size,
 323                                 in_rate,
 324                                 out_rate);
 325                 }
 326         }
 327
 328
 329         if(direction == PLAY_FORWARD)
 330                 this->output_position = out_position + out_len;
 331         else
 332                 this->output_position = out_position - out_len;
 333
 334 //printf("Resample::resample 2 %d %d\n", this->output_position, out_position);
 335 //printf("Resample::resample 2 %d %d\n", out_len, output_size);
 336
 337 //printf("Resample::resample 2 %d\n", output_size);
 338         return result;
 339 }
 340
 341
 342