14b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian/* 24cb4f7ce9184de9a221239c28afcf912e7e1ed43Dan Bornstein * Copyright (C) 2007 The Android Open Source Project 34b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * 44b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * Licensed under the Apache License, Version 2.0 (the "License"); 54b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * you may not use this file except in compliance with the License. 64b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * You may obtain a copy of the License at 74b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * 84b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * http://www.apache.org/licenses/LICENSE-2.0 94b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * 104b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * Unless required by applicable law or agreed to in writing, software 114b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * distributed under the License is distributed on an "AS IS" BASIS, 124b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 134b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * See the License for the specific language governing permissions and 144b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian * limitations under the License. 154b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian */ 164b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian 174b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian#include <math.h> 184b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian#include <stdio.h> 1973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger#include <unistd.h> 2073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger#include <stdlib.h> 2173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger#include <string.h> 224b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian 234b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopianstatic double sinc(double x) { 244b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian if (fabs(x) == 0.0f) return 1.0f; 254b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian return sin(x) / x; 264b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian} 274b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian 284b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopianstatic double sqr(double x) { 294b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian return x*x; 304b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian} 314b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian 324b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopianstatic double I0(double x) { 334b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian // from the Numerical Recipes in C p. 237 344b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian double ax,ans,y; 354b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian ax=fabs(x); 364b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian if (ax < 3.75) { 374b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian y=x/3.75; 384b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian y*=y; 394b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian ans=1.0+y*(3.5156229+y*(3.0899424+y*(1.2067492 4073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger +y*(0.2659732+y*(0.360768e-1+y*0.45813e-2))))); 414b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian } else { 424b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian y=3.75/ax; 434b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian ans=(exp(ax)/sqrt(ax))*(0.39894228+y*(0.1328592e-1 4473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger +y*(0.225319e-2+y*(-0.157565e-2+y*(0.916281e-2 4573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger +y*(-0.2057706e-1+y*(0.2635537e-1+y*(-0.1647633e-1 4673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger +y*0.392377e-2)))))))); 474b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian } 484b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian return ans; 494b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian} 504b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian 5173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflingerstatic double kaiser(int k, int N, double beta) { 524b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian if (k < 0 || k > N) 534b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian return 0; 5473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger return I0(beta * sqrt(1.0 - sqr((2.0*k)/N - 1.0))) / I0(beta); 5573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger} 5673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 5773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 5873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflingerstatic void usage(char* name) { 5973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger fprintf(stderr, 6073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger "usage: %s [-h] [-d] [-s sample_rate] [-c cut-off_frequency] [-n half_zero_crossings] [-f {float|fixed}] [-b beta] [-v dBFS] [-l lerp]\n" 6173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " %s [-h] [-d] [-s sample_rate] [-c cut-off_frequency] [-n half_zero_crossings] [-f {float|fixed}] [-b beta] [-v dBFS] -p M/N\n" 6273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -h this help message\n" 6373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -d debug, print comma-separated coefficient table\n" 6473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -p generate poly-phase filter coefficients, with sample increment M/N\n" 6573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -s sample rate (48000)\n" 6673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -c cut-off frequency (20478)\n" 6773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -n number of zero-crossings on one side (8)\n" 6873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -l number of lerping bits (4)\n" 6973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -f output format, can be fixed-point or floating-point (fixed)\n" 7073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -b kaiser window parameter beta (7.865 [-80dB])\n" 7173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger " -v attenuation in dBFS (0)\n", 7273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger name, name 7373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger ); 7473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger exit(0); 754b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian} 764b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian 774b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopianint main(int argc, char** argv) 784b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian{ 794b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian // nc is the number of bits to store the coefficients 8073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger const int nc = 32; 814b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian 8273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger bool polyphase = false; 8373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger unsigned int polyM = 160; 8473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger unsigned int polyN = 147; 8573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger bool debug = false; 8673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger double Fs = 48000; 87b4b75b47c2a4248e60bbc3229d6acc4d5f872431Mathias Agopian double Fc = 20478; 8873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger double atten = 1; 8973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger int format = 0; 904b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian 912a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian 9273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // in order to keep the errors associated with the linear 9373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // interpolation of the coefficients below the quantization error 9473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // we must satisfy: 9573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 2^nz >= 2^(nc/2) 9673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 9773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // for 16 bit coefficients that would be 256 9873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 9973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // note that increasing nz only increases memory requirements, 10073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // but doesn't increase the amount of computation to do. 10173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 10273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 10373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // see: 10473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // Smith, J.O. Digital Audio Resampling Home Page 10573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // https://ccrma.stanford.edu/~jos/resample/, 2011-03-29 10673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 10773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger int nz = 4; 10873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 10973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // | 0.1102*(A - 8.7) A > 50 11073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // beta = | 0.5842*(A - 21)^0.4 + 0.07886*(A - 21) 21 <= A <= 50 11173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // | 0 A < 21 11273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // with A is the desired stop-band attenuation in dBFS 11373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 11473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // for eg: 11573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 1162a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian // 30 dB 2.210 1172a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian // 40 dB 3.384 1182a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian // 50 dB 4.538 1192a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian // 60 dB 5.658 1202a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian // 70 dB 6.764 1212a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian // 80 dB 7.865 1222a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian // 90 dB 8.960 1232a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian // 100 dB 10.056 12473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger double beta = 7.865; 12573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 12673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 12773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 2*nzc = (A - 8) / (2.285 * dw) 12873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // with dw the transition width = 2*pi*dF/Fs 12973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 13073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger int nzc = 8; 13173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 13273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 13373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // Example: 13473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 44.1 KHz to 48 KHz resampling 13573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 100 dB rejection above 28 KHz 13673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // (the spectrum will fold around 24 KHz and we want 100 dB rejection 13773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // at the point where the folding reaches 20 KHz) 13873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // ...___|_____ 13973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // | \| 14073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // | ____/|\____ 14173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // |/alias| \ 14273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // ------/------+------\---------> KHz 14373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 20 24 28 14473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 14573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // Transition band 8 KHz, or dw = 1.0472 14673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 14773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // beta = 10.056 14873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // nzc = 20 14973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // 15073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 15173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger int ch; 15273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger while ((ch = getopt(argc, argv, ":hds:c:n:f:l:b:p:v:")) != -1) { 15373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger switch (ch) { 15473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 'd': 15573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger debug = true; 15673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 15773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 'p': 15873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (sscanf(optarg, "%u/%u", &polyM, &polyN) != 2) { 15973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger usage(argv[0]); 16073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 16173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger polyphase = true; 16273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 16373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 's': 16473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger Fs = atof(optarg); 16573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 16673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 'c': 16773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger Fc = atof(optarg); 16873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 16973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 'n': 17073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger nzc = atoi(optarg); 17173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 17273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 'l': 17373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger nz = atoi(optarg); 17473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 17573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 'f': 17673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (!strcmp(optarg,"fixed")) format = 0; 17773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger else if (!strcmp(optarg,"float")) format = 1; 17873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger else usage(argv[0]); 17973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 18073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 'b': 18173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger beta = atof(optarg); 18273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 18373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 'v': 18473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger atten = pow(10, -fabs(atof(optarg))*0.05 ); 18573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 18673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger case 'h': 18773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger default: 18873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger usage(argv[0]); 18973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger break; 19073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 19173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 19273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 19373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // cut off frequency ratio Fc/Fs 19473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger double Fcr = Fc / Fs; 19573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 1962a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian 19773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // total number of coefficients (one side) 19846afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian const int M = (1 << nz); 19946afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian const int N = M * nzc; 2004b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian 2014b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian // generate the right half of the filter 20273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (!debug) { 20373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("// cmd-line: "); 20473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger for (int i=1 ; i<argc ; i++) { 20573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("%s ", argv[i]); 20673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 20773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("\n"); 20873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (!polyphase) { 20973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("const int32_t RESAMPLE_FIR_SIZE = %d;\n", N); 21073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("const int32_t RESAMPLE_FIR_LERP_INT_BITS = %d;\n", nz); 21173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("const int32_t RESAMPLE_FIR_NUM_COEF = %d;\n", nzc); 21273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } else { 21373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("const int32_t RESAMPLE_FIR_SIZE = %d;\n", 2*nzc*polyN); 21473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("const int32_t RESAMPLE_FIR_NUM_COEF = %d;\n", 2*nzc); 21573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 21673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (!format) { 21773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("const int32_t RESAMPLE_FIR_COEF_BITS = %d;\n", nc); 21873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 21973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("\n"); 22073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("static %s resampleFIR[] = {", !format ? "int32_t" : "float"); 22173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 22273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 22373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (!polyphase) { 22446afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian for (int i=0 ; i<=M ; i++) { // an extra set of coefs for interpolation 22546afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian for (int j=0 ; j<nzc ; j++) { 22646afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian int ix = j*M + i; 22746afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian double x = (2.0 * M_PI * ix * Fcr) / (1 << nz); 22846afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian double y = kaiser(ix+N, 2*N, beta) * sinc(x) * 2.0 * Fcr; 22946afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian y *= atten; 23073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 23146afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian if (!debug) { 23246afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian if (j == 0) 23346afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian printf("\n "); 23446afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian } 23573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 23646afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian if (!format) { 23746afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian int64_t yi = floor(y * ((1ULL<<(nc-1))) + 0.5); 23846afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian if (yi >= (1LL<<(nc-1))) yi = (1LL<<(nc-1))-1; 23946afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian printf("0x%08x, ", int32_t(yi)); 24046afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian } else { 24146afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian printf("%.9g%s ", y, debug ? "," : "f,"); 24246afbec3743f1d799f185273ff897d1f8e0175ddMathias Agopian } 24373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 24473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 24573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } else { 24673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger for (int j=0 ; j<polyN ; j++) { 24773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // calculate the phase 24873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger double p = ((polyM*j) % polyN) / double(polyN); 24973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (!debug) printf("\n "); 25073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger else printf("\n"); 25173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger // generate a FIR per phase 25273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger for (int i=-nzc ; i<nzc ; i++) { 25373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger double x = 2.0 * M_PI * Fcr * (i + p); 254d88a051aff15fdf5c57e1e5a4083bbd9635af3adMathias Agopian double y = kaiser(i+N, 2*N, beta) * sinc(x) * 2.0 * Fcr;; 25573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger y *= atten; 25673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (!format) { 25773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger int64_t yi = floor(y * ((1ULL<<(nc-1))) + 0.5); 25873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (yi >= (1LL<<(nc-1))) yi = (1LL<<(nc-1))-1; 25973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("0x%08x", int32_t(yi)); 26073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } else { 26173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("%.9g%s", y, debug ? "" : "f"); 26273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 26373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 26473e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (debug && (i==nzc-1)) { 26573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } else { 26673e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf(", "); 26773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 26873e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 26973e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 27073e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger } 27173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 27273e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger if (!debug) { 27373e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("\n};"); 2744b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian } 27573e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger printf("\n"); 2764b8a3d8a89814dc3fb365f18d01733e26eb495a1Mathias Agopian return 0; 27773e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger} 2782a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian 2792a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian// http://www.csee.umbc.edu/help/sound/AFsp-V2R1/html/audio/ResampAudio.html 2802a967b3fff07b8711aef41f839ad7521323bb64dMathias Agopian 28173e90268adf4c9638b8d820a802e5e9a8ebe6597Pixelflinger 282