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