Main Page | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Class Members | File Members

resample.C

Go to the documentation of this file.
00001 #include "clip.h"
00002 #include "file.h"
00003 #include "resample.h"
00004 
00005 #include <math.h>
00006 #include <stdio.h>
00007 #include <string.h>
00008 
00009 // Resampling from Lame
00010 
00011 Resample::Resample(File *file, int channels)
00012 {
00013 //printf("Resample::Resample 1 %d\n", channels);
00014         this->file = file;
00015         this->channels = channels;
00016 
00017         old = new double*[channels];
00018         for(int i = 0; i < channels; i++)
00019         {
00020                 old[i] = new double[BLACKSIZE];
00021         }
00022         itime = new double[channels];
00023         output_temp_start = new long[channels];
00024         bzero(output_temp_start, sizeof(long) * channels);
00025         resample_init = new int[channels];
00026         bzero(resample_init, sizeof(int) * channels);
00027         last_ratio = 0;
00028         output_temp = 0;
00029         output_size = new long[channels];
00030         bzero(output_size, sizeof(long) * channels);
00031         output_allocation = 0;
00032         input_size = RESAMPLE_CHUNKSIZE;
00033         input_chunk_end = new long[channels];
00034         bzero(input_chunk_end, sizeof(long) * channels);
00035         input = new double[input_size];
00036         last_out_end = new long[channels];
00037         bzero(last_out_end, sizeof(long) * channels);
00038 //printf("Resample::Resample 2 %d\n", channels);
00039 }
00040 
00041 
00042 Resample::~Resample()
00043 {
00044         for(int i = 0; i < channels; i++)
00045         {
00046                 delete [] old[i];
00047         }
00048         if(output_temp)
00049         {
00050                 for(int i = 0; i < channels; i++)
00051                 {
00052                         delete [] output_temp[i];
00053                 }
00054                 delete [] output_temp;
00055         }
00056 
00057         delete [] input_chunk_end;
00058         delete [] input;
00059         delete [] old;
00060         delete [] itime;
00061         delete [] output_temp_start;
00062         delete [] output_size;
00063         delete [] last_out_end;
00064 }
00065 
00066 void Resample::reset(int channel)
00067 {
00068 //printf("Resample::reset 1 channel=%d normalized_sample_rate=%d\n", channel, file->normalized_sample_rate);
00069         if(channel < 0)
00070         {
00071                 bzero(resample_init, sizeof(int) * channels);
00072                 bzero(output_size, sizeof(long) * channels);
00073                 bzero(last_out_end, sizeof(long) * channels);
00074                 bzero(input_chunk_end, sizeof(long) * channels);
00075         }
00076         else
00077         {
00078                 resample_init[channel] = 0;
00079                 output_size[channel] = 0;
00080                 last_out_end[channel] = 0;
00081                 input_chunk_end[channel] = 0;
00082         }
00083 }
00084 
00085 double Resample::blackman(int i, double offset, double fcn, int l)
00086 {
00087   /* This algorithm from:
00088 SIGNAL PROCESSING ALGORITHMS IN FORTRAN AND C
00089 S.D. Stearns and R.A. David, Prentice-Hall, 1992
00090   */
00091 
00092         double bkwn;
00093         double wcn = (M_PI * fcn);
00094         double dly = l / 2.0;
00095         double x = i-offset;
00096         if(x < 0) x = 0;
00097         else
00098         if(x > l) x = l;
00099 
00100         bkwn = 0.42 - 0.5 * cos((x * 2) * M_PI /l) + 0.08 * cos((x * 4) * M_PI /l);
00101         if(fabs(x - dly) < 1e-9) 
00102                 return wcn / M_PI;
00103     else 
00104         return (sin((wcn * (x - dly))) / (M_PI * (x - dly)) * bkwn);
00105 }
00106 
00107 
00108 int Resample::get_output_size(int channel)
00109 {
00110         return output_size[channel];
00111 }
00112 
00113 void Resample::read_output(double *output, int channel, int size)
00114 {
00115         memcpy(output, output_temp[channel], size * sizeof(double));
00116 // Shift leftover forward
00117         for(int i = size; i < output_size[channel]; i++)
00118                 output_temp[channel][i - size] = output_temp[channel][i];
00119         output_size[channel] -= size;
00120 }
00121 
00122 
00123 
00124 void Resample::resample_chunk(double *input,
00125         long in_len,
00126         int in_rate,
00127         int out_rate,
00128         int channel)
00129 {
00130         double resample_ratio = (double)in_rate / out_rate;
00131         int filter_l;
00132         double fcn, intratio;
00133         double offset, xvalue;
00134         int num_used;
00135         int i, j, k;
00136 
00137         intratio = (fabs(resample_ratio - floor(.5 + resample_ratio)) < .0001);
00138         fcn = .90 / resample_ratio;
00139         if(fcn > .90) fcn = .90;
00140         filter_l = BLACKSIZE - 6;  
00141 /* must be odd */
00142         if(0 == filter_l % 2 ) --filter_l;  
00143 
00144 /* if resample_ratio = int, filter_l should be even */
00145         filter_l += (int)intratio;
00146 
00147 // Blackman filter initialization must be called whenever there is a 
00148 // sampling ratio change
00149         if(!resample_init[channel] || last_ratio != resample_ratio)
00150         {
00151                 resample_init[channel] = 1;
00152                 itime[channel] = 0;
00153                 bzero(old[channel], sizeof(double) * BLACKSIZE);
00154 
00155 // precompute blackman filter coefficients
00156         for (j = 0; j <= 2 * BPC; ++j) 
00157                 {
00158                         for(j = 0; j <= 2 * BPC; j++)
00159                         {
00160                                 offset = (double)(j - BPC) / (2 * BPC);
00161                                 for(i = 0; i <= filter_l; i++)
00162                                 {
00163                                         blackfilt[j][i] = blackman(i, offset, fcn, filter_l);
00164                                 }
00165                         }
00166                 }
00167         }
00168 
00169 // Main loop
00170         double *inbuf_old = old[channel];
00171         for(k = 0; 1; k++)
00172         {
00173                 double time0;
00174                 int joff;
00175                 
00176                 time0 = k * resample_ratio;
00177                 j = (int)floor(time0 - itime[channel]);
00178 
00179 //              if(j + filter_l / 2 >= input_size) break;
00180                 if(j + (filter_l + 1) / 2 >= in_len) break;
00181 
00182 /* blackman filter.  by default, window centered at j+.5(filter_l%2) */
00183 /* but we want a window centered at time0.   */
00184                 offset = (time0 - itime[channel] - (j + .5 * (filter_l % 2)));
00185                 joff = (int)floor((offset * 2 * BPC) + BPC + .5);
00186                 xvalue = 0;
00187 
00188                 for(i = 0; i <= filter_l; i++)
00189                 {
00190                         int j2 = i + j - filter_l / 2;
00191                         double y = ((j2 < 0) ? inbuf_old[BLACKSIZE + j2] : input[j2]);
00192 
00193                         xvalue += y * blackfilt[joff][i];
00194                 }
00195                 
00196                 if(output_allocation <= output_size[channel])
00197                 {
00198                         double **new_output = new double*[channels];
00199                         long new_allocation = output_allocation ? (output_allocation * 2) : 16384;
00200                         for(int l = 0; l < channels; l++)
00201                         {
00202                                 new_output[l] = new double[new_allocation];
00203                                 if(output_temp) 
00204                                 {
00205                                         bcopy(output_temp[l], new_output[l], output_allocation * sizeof(double));
00206                                         delete [] output_temp[l];
00207                                 }
00208                         }
00209 
00210                         if(output_temp) delete [] output_temp;
00211                         output_temp = new_output;
00212                         output_allocation = new_allocation;
00213                 }
00214 
00215                 output_temp[channel][output_size[channel]++] = xvalue;
00216         }
00217 
00218         num_used = MIN(in_len, j + filter_l / 2);
00219         itime[channel] += num_used - k * resample_ratio;
00220 //      for(i = 0; i < BLACKSIZE; i++)
00221 //              inbuf_old[i] = input[num_used + i - BLACKSIZE];
00222         bcopy(input + num_used - BLACKSIZE, inbuf_old, BLACKSIZE * sizeof(double));
00223 
00224         last_ratio = resample_ratio;
00225 }
00226 
00227 void Resample::read_chunk(double *input, long len, int &reseek, int iteration)
00228 {
00229 //printf("Resample::read_chunk 1\n");
00230         if(reseek)
00231         {
00232                 file->set_audio_position(file->current_sample, 0);
00233                 reseek= 0;
00234         }
00235         else
00236         if(iteration == 0)
00237         {
00238 // Resume at the end of the last resample call
00239                 file->set_audio_position(input_chunk_end[file->current_channel], 0);
00240         }
00241 
00242         file->read_samples(input, len, 0);
00243         input_chunk_end[file->current_channel] = file->current_sample;
00244 
00245 //printf("Resample::read_chunk 2\n");
00246 }
00247 
00248 int Resample::resample(double *output, 
00249         long out_len,
00250         int in_rate,
00251         int out_rate,
00252         int channel,
00253         long in_position,
00254         long out_position)
00255 {
00256         int total_input = 0;
00257         int reseek = 0;
00258 
00259 #define REPOSITION(x, y) \
00260         (labs((x) - (y)) > 1)
00261 
00262 
00263 
00264 //printf("Resample::resample 1 last_out_end=%d out_position=%d\n", last_out_end[channel], out_position);
00265 
00266         if(REPOSITION(last_out_end[channel], out_position))
00267         {
00268                 reseek = 1;
00269                 reset(channel);
00270         }
00271 
00272 
00273 
00274 
00275 
00276 
00277         output_temp_start[channel] = file->get_audio_position(out_rate) + out_len;
00278         last_out_end[channel] = out_position + out_len;
00279 
00280         int i = 0;
00281         while(out_len > 0)
00282         {
00283 // Drain output buffer
00284                 if(output_size[channel])
00285                 {
00286                         int fragment_len = output_size[channel];
00287                         if(fragment_len > out_len) fragment_len = out_len;
00288 
00289 //printf("Resample::resample 1 %d %d %d\n", out_len, output_size[channel], channel);
00290                         bcopy(output_temp[channel], output, fragment_len * sizeof(double));
00291 
00292 // Shift leftover forward
00293 //                      for(int i = fragment_len; i < output_size[channel]; i++)
00294 //                              output_temp[channel][i - fragment_len] = output_temp[channel][i];
00295                         bcopy(output_temp[channel] + fragment_len, output_temp[channel], (output_size[channel] - fragment_len) * sizeof(double)); 
00296 
00297                         output_size[channel] -= fragment_len;
00298                         out_len -= fragment_len;
00299                         output += fragment_len;
00300                 }
00301 
00302 // Import new samples
00303 //printf("Resample::resample 2 %d %d\n", out_len, channel);
00304                 if(out_len > 0)
00305                 {
00306 //printf("Resample::resample 3 input_size=%d reseek=%d out_position=%d channel=%d\n", input_size, reseek, out_position, channel);
00307                         read_chunk(input, input_size, reseek, i);
00308                         resample_chunk(input,
00309                                 input_size,
00310                                 in_rate,
00311                                 out_rate,
00312                                 channel);
00313                         total_input += input_size;
00314                 }
00315 
00316                 i++;
00317         }
00318 //printf("Resample::resample 2 %d %d\n", last_out_end[channel], out_position);
00319 //printf("Resample::resample 2 %d %d %d\n", out_len, output_size[channel], channel);
00320 
00321 //printf("Resample::resample 2 %d %d\n", channel, output_size[channel]);
00322 
00323         return total_input;
00324 }
00325 
00326 
00327 Resample_float::Resample_float(File *file, int channels)
00328 {
00329 //printf("Resample_float::Resample_float 1 %d\n", channels);
00330         this->file = file;
00331         this->channels = channels;
00332 
00333         old = new float*[channels];
00334         for(int i = 0; i < channels; i++)
00335         {
00336                 old[i] = new float[BLACKSIZE];
00337         }
00338         itime = new float[channels];
00339         output_temp_start = new long[channels];
00340         bzero(output_temp_start, sizeof(long) * channels);
00341         resample_init = new int[channels];
00342         bzero(resample_init, sizeof(int) * channels);
00343         last_ratio = 0;
00344         output_temp = 0;
00345         output_size = new long[channels];
00346         bzero(output_size, sizeof(long) * channels);
00347         output_allocation = 0;
00348         input_size = RESAMPLE_CHUNKSIZE;
00349         input_chunk_end = new long[channels];
00350         bzero(input_chunk_end, sizeof(long) * channels);
00351         input = new float[input_size];
00352         last_out_end = new long[channels];
00353         bzero(last_out_end, sizeof(long) * channels);
00354 }
00355 
00356 
00357 Resample_float::~Resample_float()
00358 {
00359         for(int i = 0; i < channels; i++)
00360         {
00361                 delete [] old[i];
00362         }
00363         if(output_temp)
00364         {
00365                 for(int i = 0; i < channels; i++)
00366                 {
00367                         delete [] output_temp[i];
00368                 }
00369                 delete [] output_temp;
00370         }
00371 
00372         delete [] input_chunk_end;
00373         delete [] input;
00374         delete [] old;
00375         delete [] itime;
00376         delete [] output_temp_start;
00377         delete [] output_size;
00378         delete [] last_out_end;
00379 }
00380 
00381 void Resample_float::reset(int channel)
00382 {
00383 //printf("Resample_float::reset 1 channel=%d normalized_sample_rate=%d\n", channel, file->normalized_sample_rate);
00384         if(channel < 0)
00385         {
00386                 bzero(resample_init, sizeof(int) * channels);
00387                 bzero(output_size, sizeof(long) * channels);
00388                 bzero(last_out_end, sizeof(long) * channels);
00389                 bzero(input_chunk_end, sizeof(long) * channels);
00390         }
00391         else
00392         {
00393                 resample_init[channel] = 0;
00394                 output_size[channel] = 0;
00395                 last_out_end[channel] = 0;
00396                 input_chunk_end[channel] = 0;
00397         }
00398 }
00399 
00400 float Resample_float::blackman(int i, float offset, float fcn, int l)
00401 {
00402   /* This algorithm from:
00403 SIGNAL PROCESSING ALGORITHMS IN FORTRAN AND C
00404 S.D. Stearns and R.A. David, Prentice-Hall, 1992
00405   */
00406 
00407         float bkwn;
00408         float wcn = (M_PI * fcn);
00409         float dly = l / 2.0;
00410         float x = i-offset;
00411         if(x < 0) x = 0;
00412         else
00413         if(x > l) x = l;
00414 
00415         bkwn = 0.42 - 0.5 * cos((x * 2) * M_PI /l) + 0.08 * cos((x * 4) * M_PI /l);
00416         if(fabs(x - dly) < 1e-9) 
00417                 return wcn / M_PI;
00418     else 
00419         return (sin((wcn * (x - dly))) / (M_PI * (x - dly)) * bkwn);
00420 }
00421 
00422 
00423 int Resample_float::get_output_size(int channel)
00424 {
00425         return output_size[channel];
00426 }
00427 
00428 void Resample_float::read_output(double *output, int channel, int size)
00429 {
00430         memcpy(output, output_temp[channel], size * sizeof(double));
00431 // Shift leftover forward
00432         for(int i = size; i < output_size[channel]; i++)
00433                 output_temp[channel][i - size] = output_temp[channel][i];
00434         output_size[channel] -= size;
00435 }
00436 
00437 
00438 
00439 void Resample_float::resample_chunk(float *input,
00440         long in_len,
00441         int in_rate,
00442         int out_rate,
00443         int channel)
00444 {
00445         float resample_ratio = (float)in_rate / out_rate;
00446         int filter_l;
00447         float fcn, intratio;
00448         float offset, xvalue;
00449         int num_used;
00450         int i, j, k;
00451 
00452 //printf("Resample_float::resample_chunk 1\n");
00453         intratio = (fabs(resample_ratio - floor(.5 + resample_ratio)) < .0001);
00454         fcn = .90 / resample_ratio;
00455         if(fcn > .90) fcn = .90;
00456         filter_l = BLACKSIZE - 6;  
00457 /* must be odd */
00458         if(0 == filter_l % 2 ) --filter_l;  
00459 
00460 //printf("Resample_float::resample_chunk 2\n");
00461 /* if resample_ratio = int, filter_l should be even */
00462         filter_l += (int)intratio;
00463 
00464 //printf("Resample_float::resample_chunk 3\n");
00465 // Blackman filter initialization must be called whenever there is a 
00466 // sampling ratio change
00467         if(!resample_init[channel] || last_ratio != resample_ratio)
00468         {
00469                 resample_init[channel] = 1;
00470                 itime[channel] = 0;
00471                 bzero(old[channel], sizeof(float) * BLACKSIZE);
00472 
00473 //printf("Resample_float::resample_chunk 4\n");
00474 // precompute blackman filter coefficients
00475         for (j = 0; j <= 2 * BPC; ++j) 
00476                 {
00477                         for(j = 0; j <= 2 * BPC; j++)
00478                         {
00479                                 offset = (float)(j - BPC) / (2 * BPC);
00480                                 for(i = 0; i <= filter_l; i++)
00481                                 {
00482                                         blackfilt[j][i] = blackman(i, offset, fcn, filter_l);
00483                                 }
00484                         }
00485                 }
00486         }
00487 
00488 //printf("Resample_float::resample_chunk 5\n");
00489 // Main loop
00490         float *inbuf_old = old[channel];
00491         for(k = 0; 1; k++)
00492         {
00493                 float time0;
00494                 int joff;
00495                 
00496 //printf("Resample_float::resample_chunk 6\n");
00497                 time0 = k * resample_ratio;
00498                 j = (int)floor(time0 - itime[channel]);
00499 
00500 //              if(j + filter_l / 2 >= input_size) break;
00501                 if(j + (filter_l + 1) / 2 >= in_len) break;
00502 
00503 //printf("Resample_float::resample_chunk 7\n");
00504 /* blackman filter.  by default, window centered at j+.5(filter_l%2) */
00505 /* but we want a window centered at time0.   */
00506                 offset = (time0 - itime[channel] - (j + .5 * (filter_l % 2)));
00507                 joff = (int)floor((offset * 2 * BPC) + BPC + .5);
00508                 xvalue = 0;
00509 
00510 //printf("Resample_float::resample_chunk 8\n");
00511                 for(i = 0; i <= filter_l; i++)
00512                 {
00513                         int j2 = i + j - filter_l / 2;
00514                         float y = ((j2 < 0) ? inbuf_old[BLACKSIZE + j2] : input[j2]);
00515 
00516 //printf("Resample_float::resample_chunk 9\n");
00517                         xvalue += (double)y * blackfilt[joff][i];
00518                 }
00519                 
00520 //printf("Resample_float::resample_chunk 10\n");
00521                 if(output_allocation <= output_size[channel])
00522                 {
00523                         double **new_output = new double*[channels];
00524                         long new_allocation = output_allocation ? (output_allocation * 2) : 16384;
00525                         for(int l = 0; l < channels; l++)
00526                         {
00527                                 new_output[l] = new double[new_allocation];
00528                                 if(output_temp) 
00529                                 {
00530                                         bcopy(output_temp[l], new_output[l], output_allocation * sizeof(double));
00531                                         delete [] output_temp[l];
00532                                 }
00533                         }
00534 
00535                         if(output_temp) delete [] output_temp;
00536                         output_temp = new_output;
00537                         output_allocation = new_allocation;
00538                 }
00539 
00540 //printf("Resample_float::resample_chunk 11 %d %d\n", output_size[channel], output_allocation);
00541                 output_temp[channel][output_size[channel]++] = xvalue;
00542         }
00543 
00544 //printf("Resample_float::resample_chunk 12\n");
00545         num_used = MIN(in_len, j + filter_l / 2);
00546         itime[channel] += num_used - k * resample_ratio;
00547 //      for(i = 0; i < BLACKSIZE; i++)
00548 //              inbuf_old[i] = input[num_used + i - BLACKSIZE];
00549         bcopy(input + num_used - BLACKSIZE, inbuf_old, BLACKSIZE * sizeof(float));
00550 
00551 //printf("Resample_float::resample_chunk 13\n");
00552         last_ratio = resample_ratio;
00553 }
00554 
00555 void Resample_float::read_chunk(float *input, long len, int &reseek, int iteration)
00556 {
00557 //printf("Resample_float::read_chunk 1\n");
00558         if(reseek)
00559         {
00560                 file->set_audio_position(file->current_sample, 0);
00561                 reseek= 0;
00562         }
00563         else
00564         if(iteration == 0)
00565         {
00566 // Resume at the end of the last resample call
00567                 file->set_audio_position(input_chunk_end[file->current_channel], 0);
00568         }
00569 
00570         file->read_samples(0, len, 0, input);
00571         input_chunk_end[file->current_channel] = file->current_sample;
00572 
00573 //printf("Resample_float::read_chunk 2\n");
00574 }
00575 
00576 int Resample_float::resample(double *output, 
00577         long out_len,
00578         int in_rate,
00579         int out_rate,
00580         int channel,
00581         long in_position,
00582         long out_position)
00583 {
00584         int total_input = 0;
00585         int reseek = 0;
00586 
00587 #define REPOSITION(x, y) \
00588         (labs((x) - (y)) > 1)
00589 
00590 
00591 
00592 //printf("Resample_float::resample 1 last_out_end=%d out_position=%d\n", last_out_end[channel], out_position);
00593 
00594         if(REPOSITION(last_out_end[channel], out_position))
00595         {
00596                 reseek = 1;
00597                 reset(channel);
00598         }
00599 
00600         output_temp_start[channel] = file->get_audio_position(out_rate) + out_len;
00601         last_out_end[channel] = out_position + out_len;
00602 
00603         int i = 0;
00604         while(out_len > 0)
00605         {
00606 // Drain output buffer
00607                 if(output_size[channel])
00608                 {
00609                         int fragment_len = output_size[channel];
00610                         if(fragment_len > out_len) fragment_len = out_len;
00611 
00612 //printf("Resample_float::resample 1 %d %d %d\n", out_len, output_size[channel], channel);
00613                         bcopy(output_temp[channel], output, fragment_len * sizeof(double));
00614 
00615 
00616 // Shift leftover forward
00617                         //for(int i = fragment_len; i < output_size[channel]; i++)
00618                         //      output_temp[channel][i - fragment_len] = output_temp[channel][i];
00619                         bcopy(output_temp[channel] + fragment_len, output_temp[channel], (output_size[channel] - fragment_len) * sizeof(double)); 
00620 
00621                         output_size[channel] -= fragment_len;
00622                         out_len -= fragment_len;
00623                         output += fragment_len;
00624                 }
00625 
00626 // Import new samples
00627 //printf("Resample_float::resample 2 %d %d\n", out_len, channel);
00628                 if(out_len > 0)
00629                 {
00630 //printf("Resample_float::resample 3 input_size=%d reseek=%d out_position=%d channel=%d\n", input_size, reseek, out_position, channel);
00631                         read_chunk(input, input_size, reseek, i);
00632                         resample_chunk(input,
00633                                 input_size,
00634                                 in_rate,
00635                                 out_rate,
00636                                 channel);
00637                         total_input += input_size;
00638                 }
00639 
00640                 i++;
00641         }
00642 //printf("Resample_float::resample 2 %d %d\n", last_out_end[channel], out_position);
00643 //printf("Resample_float::resample 2 %d %d %d\n", out_len, output_size[channel], channel);
00644 
00645 //printf("Resample_float::resample 2 %d %d\n", channel, output_size[channel]);
00646 
00647         return total_input;
00648 }

Generated on Sun Jan 8 13:39:00 2006 for Cinelerra-svn by  doxygen 1.4.4