hvirtual/cinelerra/resample.C

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