窗函数的C语言实现
来源：互联网 发布：阿里云cdn防doss 编辑：程序博客网 时间：2024/06/05 22:36
一般的讲数字信号处理的书中都会提到窗函数。大多数只会提及其中的几种。这里我把这些窗都用C语言实现了一下，都不复杂，但如果要自己去弄也挺费时间。所有函数都用Matlab验证了。包括以下窗：
WindowFunction.hWindowFunction.cWindowFunction.h
1 /*窗类型*/ 2 typedef enum 3 { 4 Bartlett = 0, 5 BartLettHann, 6 BlackMan, 7 BlackManHarris, 8 Bohman, 9 Chebyshev,10 FlatTop,11 Gaussian,12 Hamming,13 Hann,14 Kaiser,15 Nuttal,16 Parzen,17 Rectangular,18 Taylor, 19 Triangular,20 Tukey21 }winType;
别的不多说了，直接上干货。
1 /* 2 *file WindowFunction.h 3 *author Vincent Cui 4 *e-mail whcui1987@163.com 5 *version 0.3 6 *data 31-Oct-2014 7 *brief 各种窗函数的C语言实现 8 */ 9 10 11 #ifndef _WINDOWFUNCTION_H_12 #define _WINDOWFUNCTION_H_13 14 #include "GeneralConfig.h"15 16 #define BESSELI_K_LENGTH 1017 18 #define FLATTOPWIN_A0 0.21557899519 #define FLATTOPWIN_A1 0.4166315820 #define FLATTOPWIN_A2 0.27726315821 #define FLATTOPWIN_A3 0.08357894722 #define FLATTOPWIN_A4 0.00694736823 24 #define NUTTALL_A0 0.363581925 #define NUTTALL_A1 0.489177526 #define NUTTALL_A2 0.136599527 #define NUTTALL_A3 0.010641128 29 #define BLACKMANHARRIS_A0 0.3587530 #define BLACKMANHARRIS_A1 0.4882931 #define BLACKMANHARRIS_A2 0.1412832 #define BLACKMANHARRIS_A3 0.0116833 34 35 36 dspErrorStatus taylorWin(dspUint_16 N, dspUint_16 nbar, dspDouble sll, dspDouble **w);37 dspErrorStatus triangularWin(dspUint_16 N, dspDouble **w);38 dspErrorStatus tukeyWin(dspUint_16 N, dspDouble r, dspDouble **w);39 dspErrorStatus bartlettWin(dspUint_16 N, dspDouble **w);40 dspErrorStatus bartLettHannWin(dspUint_16 N, dspDouble **w);41 dspErrorStatus blackManWin(dspUint_16 N, dspDouble **w);42 dspErrorStatus blackManHarrisWin(dspUint_16 N, dspDouble **w);43 dspErrorStatus bohmanWin(dspUint_16 N, dspDouble **w);44 dspErrorStatus chebyshevWin(dspUint_16 N, dspDouble r, dspDouble **w);45 dspErrorStatus flatTopWin(dspUint_16 N, dspDouble **w);46 dspErrorStatus gaussianWin(dspUint_16 N, dspDouble alpha, dspDouble **w);47 dspErrorStatus hammingWin(dspUint_16 N, dspDouble **w);48 dspErrorStatus hannWin(dspUint_16 N, dspDouble **w);49 dspErrorStatus kaiserWin(dspUint_16 N, dspDouble beta, dspDouble **w);50 dspErrorStatus nuttalWin(dspUint_16 N, dspDouble **w);51 dspErrorStatus parzenWin(dspUint_16 N, dspDouble **w);52 dspErrorStatus rectangularWin(dspUint_16 N, dspDouble **w);53 54 55 56 #endif
1 /* 2 *file WindowFunction.c 3 *author Vincent Cui 4 *e-mail whcui1987@163.com 5 *version 0.3 6 *data 31-Oct-2014 7 *brief 各种窗函数的C语言实现 8 */ 9 10 11 #include "WindowFunction.h" 12 #include "GeneralConfig.h" 13 #include "MathReplenish.h" 14 #include "math.h" 15 #include <stdlib.h> 16 #include <stdio.h> 17 #include <string.h> 18 19 /*函数名：taylorWin 20 *说明：计算泰勒窗。泰勒加权函数 21 *输入： 22 *输出： 23 *返回： 24 *调用：prod()连乘函数 25 *其它：用过以后，需要手动释放掉*w的内存空间 26 * 调用示例：ret = taylorWin(99, 4, 40, &w); 注意此处的40是正数 表示-40dB 27 */ 28 dspErrorStatus taylorWin(dspUint_16 N, dspUint_16 nbar, dspDouble sll, dspDouble **w) 29 { 30 dspDouble A; 31 dspDouble *retDspDouble; 32 dspDouble *sf; 33 dspDouble *result; 34 dspDouble alpha,beta,theta; 35 dspUint_16 i,j; 36 37 /*A = R cosh(PI, A) = R*/ 38 A = (dspDouble)acosh(pow((dspDouble)10.0,(dspDouble)sll/20.0)) / PI; 39 A = A * A; 40 41 /*开出存放系数的空间*/ 42 retDspDouble = (dspDouble *)malloc(sizeof(dspDouble) * (nbar - 1)); 43 if(retDspDouble == NULL) 44 return DSP_ERROR; 45 sf = retDspDouble; 46 47 /*开出存放系数的空间*/ 48 retDspDouble = (dspDouble *)malloc(sizeof(dspDouble) * N); 49 if(retDspDouble == NULL) 50 return DSP_ERROR; 51 result = retDspDouble; 52 53 alpha = prod(1, 1, (nbar - 1)); 54 alpha *= alpha; 55 beta = (dspDouble)nbar / sqrt( A + pow((nbar - 0.5), 2) ); 56 for(i = 1; i <= (nbar - 1); i++) 57 { 58 *(sf + i - 1) = prod(1,1,(nbar -1 + i)) * prod(1,1,(nbar -1 - i)); 59 theta = 1; 60 for(j = 1; j <= (nbar - 1); j++) 61 { 62 theta *= 1 - (dspDouble)(i * i) / ( beta * beta * ( A + (j - 0.5) * (j - 0.5)) ); 63 } 64 *(sf + i - 1) = alpha * (dspDouble)theta / (*(sf + i - 1)); 65 } 66 67 /*奇数阶*/ 68 if((N % 2) == 1) 69 { 70 for(i = 0; i < N; i++) 71 { 72 alpha = 0; 73 for(j = 1; j <= (nbar - 1); j++) 74 { 75 alpha += (*(sf + j - 1)) * cos( 2 * PI * j * (dspDouble)(i - ((N-1)/2))/N ); 76 } 77 *(result + i) = 1 + 2 * alpha; 78 } 79 } 80 /*偶数阶*/ 81 else 82 { 83 for(i = 0; i < N; i++) 84 { 85 alpha = 0; 86 for(j = 1; j <= (nbar - 1); j++) 87 { 88 alpha += (*(sf + j - 1)) * cos( PI * j * (dspDouble)(2 * (i - (N/2)) + 1) / N ); 89 } 90 *(result + i) = 1 + 2 * alpha; 91 92 } 93 } 94 *w = result; 95 free(sf); 96 97 return DSP_SUCESS; 98 99 }100 101 102 /*103 *函数名：triangularWin104 *说明：计算三角窗函数105 *输入：106 *输出：107 *返回：108 *调用：109 *其它：用过以后，需要手动释放掉*w的内存空间110 * 调用示例：ret = triangularWin(99, &w); 111 */112 dspErrorStatus triangularWin(dspUint_16 N, dspDouble **w)113 {114 dspDouble *ptr;115 dspUint_16 i;116 117 ptr = (dspDouble *)malloc(N * sizeof(dspDouble));118 if(ptr == NULL)119 return DSP_ERROR;120 121 122 /*阶数为奇*/123 if((N % 2) == 1)124 {125 for(i = 0; i < ((N - 1)/2); i++)126 {127 *(ptr + i) = 2 * (dspDouble)(i + 1) / (N + 1);128 }129 for(i = ((N - 1)/2); i < N; i++)130 {131 *(ptr + i) = 2 * (dspDouble)(N - i) / (N + 1);132 }133 }134 /*阶数为偶*/135 else136 {137 for(i = 0; i < (N/2); i++)138 {139 *(ptr + i) = (i + i + 1) * (dspDouble)1 / N;140 }141 for(i = (N/2); i < N; i++)142 {143 *(ptr + i) = *(ptr + N - 1 - i);144 }145 }146 *w = ptr;147 148 return DSP_SUCESS;149 }150 151 /*152 *函数名：tukeyWin153 *说明：计算tukey窗函数154 *输入：155 *输出：156 *返回：linSpace()157 *调用：158 *其它：用过以后，需要手动释放掉*w的内存空间159 * 调用示例：ret = tukeyWin(99, 0.5, &w); 160 */161 dspErrorStatus tukeyWin(dspUint_16 N, dspDouble r, dspDouble **w)162 {163 dspErrorStatus retErrorStatus;164 dspUint_16 index;165 dspDouble *x,*result,*retPtr;166 dspDouble alpha;167 168 retErrorStatus = linSpace(0, 1, N, &x);169 if(retErrorStatus == DSP_ERROR)170 return DSP_ERROR;171 172 result = (dspDouble *)malloc(N * sizeof(dspDouble));173 if(result == NULL)174 return DSP_ERROR;175 176 /*r <= 0 就是矩形窗*/177 if(r <= 0)178 {179 retErrorStatus = rectangularWin(N, &retPtr);180 if(retErrorStatus == DSP_ERROR)181 return DSP_ERROR;182 /*将数据拷出来以后，释放调用的窗函数的空间*/183 memcpy(result, retPtr, ( N * sizeof(dspDouble)));184 free(retPtr);185 }186 /*r >= 1 就是汉宁窗*/187 else if(r >= 1)188 {189 retErrorStatus = hannWin(N, &retPtr);190 if(retErrorStatus == DSP_ERROR)191 return DSP_ERROR;192 /*将数据拷出来以后，释放调用的窗函数的空间*/193 memcpy(result, retPtr, ( N * sizeof(dspDouble)));194 free(retPtr);195 }196 else197 {198 for(index = 0; index < N; index++)199 {200 alpha = *(x + index);201 if(alpha < (r/2))202 {203 *(result + index) = (dspDouble)(1 + cos( 2 * PI * (dspDouble)(alpha - (dspDouble)r/2)/r))/2;204 }205 else if((alpha >= (r/2)) && (alpha <(1 - r/2)))206 {207 *(result + index) = 1;208 }209 else210 {211 *(result + index) = (dspDouble)(1 + cos( 2 * PI * (dspDouble)(alpha - 1 + (dspDouble)r/2)/r))/2;212 }213 214 }215 }216 217 free(x);218 219 *w = result;220 221 return DSP_SUCESS;222 223 }224 225 /*226 *函数名：bartlettWin227 *说明：计算bartlettWin窗函数228 *输入：229 *输出：230 *返回：231 *调用：232 *其它：用过以后，需要手动释放掉*w的内存空间233 * 调用示例：ret = bartlettWin(99, &w); 234 */235 dspErrorStatus bartlettWin(dspUint_16 N, dspDouble **w)236 {237 dspDouble *ret;238 dspUint_16 n;239 240 ret = (dspDouble *)malloc(N * sizeof(dspDouble));241 if(ret == NULL)242 return DSP_ERROR;243 244 for(n = 0; n < ( N - 1 ) / 2; n++)245 {246 *(ret + n) = 2 * (dspDouble)n / (N - 1);247 }248 249 for(n = ( N - 1 ) / 2; n < N; n++)250 {251 *(ret + n) = 2 - 2 * (dspDouble)n / (( N - 1 ));252 }253 254 *w = ret;255 256 return DSP_SUCESS;257 }258 259 /*260 *函数名：bartLettHannWin261 *说明：计算bartLettHannWin窗函数262 *输入：263 *输出：264 *返回：265 *调用：266 *其它：用过以后，需要手动释放掉*w的内存空间267 * 调用示例：ret = bartLettHannWin(99, &w); 268 */269 dspErrorStatus bartLettHannWin(dspUint_16 N, dspDouble **w)270 {271 dspUint_16 n;272 dspDouble *ret;273 274 ret = (dspDouble *)malloc(N * sizeof(dspDouble));275 if(ret == NULL)276 return DSP_ERROR;277 /*奇*/278 if(( N % 2 ) == 1)279 {280 for(n = 0; n < N; n++)281 {282 *(ret + n) = 0.62 - 0.48 * myAbs( ( (dspDouble)n / ( N - 1 ) ) - 0.5 ) + 0.38 * cos( 2 * PI * ( ((dspDouble)n / ( N - 1 ) ) - 0.5 ) );283 }284 for(n = 0; n < (N-1)/2; n++)285 {286 *(ret + n) = *(ret + N - 1 - n);287 }288 }289 /*偶*/290 else291 {292 for(n = 0; n < N; n++)293 {294 *(ret + n) = 0.62 - 0.48 * myAbs( ( (dspDouble)n / ( N - 1 ) ) - 0.5 ) + 0.38 * cos( 2 * PI * ( ((dspDouble)n / ( N - 1 ) ) - 0.5 ) );295 }296 for(n = 0; n < N/2; n++)297 {298 *(ret + n) = *(ret + N -1 - n);299 }300 }301 302 *w = ret;303 304 return DSP_SUCESS;305 306 }307 308 /*309 *函数名：blackManWin310 *说明：计算blackManWin窗函数311 *输入：312 *输出：313 *返回：314 *调用：315 *其它：用过以后，需要手动释放掉*w的内存空间316 * 调用示例：ret = blackManWin(99, &w); 317 */318 dspErrorStatus blackManWin(dspUint_16 N, dspDouble **w)319 {320 dspUint_16 n;321 dspDouble *ret;322 ret = (dspDouble *)malloc(N * sizeof(dspDouble));323 if(ret == NULL)324 return DSP_ERROR;325 326 for(n = 0; n < N; n++)327 {328 *(ret + n) = 0.42 - 0.5 * cos(2 * PI * (dspDouble)n / ( N - 1 )) + 0.08 * cos( 4 * PI * ( dspDouble )n / ( N - 1 ) );329 }330 331 *w = ret;332 333 return DSP_SUCESS;334 }335 336 /*337 *函数名：blackManHarrisWin338 *说明：计算blackManHarrisWin窗函数339 *输入：340 *输出：341 *返回：342 *调用：343 *其它：用过以后，需要手动释放掉*w的内存空间344 * 调用示例：ret = blackManHarrisWin(99, &w); 345 * minimum 4-term Blackman-harris window -- From Matlab346 */347 dspErrorStatus blackManHarrisWin(dspUint_16 N, dspDouble **w)348 {349 dspUint_16 n;350 dspDouble *ret;351 352 ret = (dspDouble *)malloc(N * sizeof(dspDouble));353 if(ret == NULL)354 return DSP_ERROR;355 356 for(n = 0; n < N; n++)357 {358 *(ret + n) = BLACKMANHARRIS_A0 - BLACKMANHARRIS_A1 * cos( 2 * PI * (dspDouble)n / (N) ) + \359 BLACKMANHARRIS_A2 * cos( 4 * PI * (dspDouble)n/ (N) ) - \360 BLACKMANHARRIS_A3 * cos( 6 * PI * (dspDouble)n/ (N) );361 }362 363 *w = ret;364 365 return DSP_SUCESS;366 }367 368 /*369 *函数名：bohmanWin370 *说明：计算bohmanWin窗函数371 *输入：372 *输出：373 *返回：374 *调用：375 *其它：用过以后，需要手动释放掉*w的内存空间376 * 调用示例：ret = bohmanWin(99, &w); 377 */378 dspErrorStatus bohmanWin(dspUint_16 N, dspDouble **w)379 {380 dspUint_16 n;381 dspDouble *ret;382 dspDouble x;383 ret = (dspDouble *)malloc(N * sizeof(dspDouble));384 if(ret == NULL)385 return DSP_ERROR;386 387 for(n = 0; n < N; n++)388 {389 x = -1 + n * (dspDouble)2 / ( N - 1 ) ;390 /*取绝对值*/391 x = x >= 0 ? x : ( x * ( -1 ) );392 *(ret + n) = ( 1 - x ) * cos( PI * x) + (dspDouble)(1 / PI) * sin( PI * x);393 }394 395 *w = ret;396 397 return DSP_SUCESS;398 }399 400 /*401 *函数名：chebyshevWin402 *说明：计算chebyshevWin窗函数403 *输入：404 *输出：405 *返回：406 *调用：407 *其它：用过以后，需要手动释放掉*w的内存空间408 * 调用示例：ret = chebyshevWin(99,100, &w); 409 */410 dspErrorStatus chebyshevWin(dspUint_16 N, dspDouble r, dspDouble **w)411 {412 dspUint_16 n,index;413 dspDouble *ret;414 dspDouble x, alpha, beta, theta, gama;415 416 ret = (dspDouble *)malloc(N * sizeof(dspDouble));417 if(ret == NULL)418 return DSP_ERROR;419 420 421 /*10^(r/20)*/422 theta = pow((dspDouble)10, (dspDouble)(myAbs(r)/20));423 beta = pow(cosh(acosh(theta)/(N - 1)),2);424 alpha = 1 - (dspDouble)1 / beta;425 426 if((N % 2) == 1)427 {428 /*计算一半的区间*/429 for( n = 1; n < ( N + 1 ) / 2; n++ )430 {431 gama = 1;432 for(index = 1; index < n; index++)433 {434 x = index * (dspDouble)( N - 1 - 2 * n + index) /(( n - index ) * (n + 1 -index));435 gama = gama * alpha * x + 1;436 }437 *(ret + n) = (N - 1) * alpha * gama; 438 }439 440 theta = *( ret + (N - 1)/2 );441 *ret = 1;442 443 for(n = 0; n < ( N + 1 ) / 2; n++ )444 {445 *(ret + n) = (dspDouble)(*(ret + n)) / theta;446 }447 448 /*填充另一半*/449 for(; n < N; n++)450 {451 *(ret + n) = ret[N - n - 1];452 }453 }454 else455 {456 /*计算一半的区间*/457 for( n = 1; n < ( N + 1 ) / 2; n++ )458 {459 gama = 1;460 for(index = 1; index < n; index++)461 {462 x = index * (dspDouble)( N - 1 - 2 * n + index) /(( n - index ) * (n + 1 -index));463 gama = gama * alpha * x + 1;464 }465 *(ret + n) = (N - 1) * alpha * gama; 466 }467 468 theta = *( ret + (N/2) - 1);469 *ret = 1;470 471 for(n = 0; n < ( N + 1 ) / 2; n++ )472 {473 *(ret + n) = (dspDouble)(*(ret + n)) / theta;474 }475 476 /*填充另一半*/477 for(; n < N; n++)478 {479 *(ret + n) = ret[N - n - 1];480 }481 }482 483 484 *w = ret;485 486 return DSP_SUCESS;487 }488 489 /*490 *函数名：flatTopWin491 *说明：计算flatTopWin窗函数492 *输入：493 *输出：494 *返回：495 *调用：496 *其它：用过以后，需要手动释放掉*w的内存空间497 * 调用示例：ret = flatTopWin(99, &w); 498 */499 dspErrorStatus flatTopWin(dspUint_16 N, dspDouble **w)500 {501 dspUint_16 n;502 dspDouble *ret;503 ret = (dspDouble *)malloc(N * sizeof(dspDouble));504 if(ret == NULL)505 return DSP_ERROR;506 507 for(n = 0; n < N; n++)508 {509 *(ret + n) = FLATTOPWIN_A0 - FLATTOPWIN_A1 * cos(2 * PI * (dspDouble)n / (N - 1)) +\510 FLATTOPWIN_A2 * cos(4 * PI * (dspDouble)n / (N - 1)) -\511 FLATTOPWIN_A3 * cos(6 * PI * (dspDouble)n / (N - 1)) +\512 FLATTOPWIN_A4 * cos(8 * PI * (dspDouble)n / (N - 1));513 }514 515 *w = ret;516 517 return DSP_SUCESS;518 }519 520 521 /*522 *函数名：gaussianWin523 *说明：计算gaussianWin窗函数524 *输入：525 *输出：526 *返回：527 *调用：528 *其它：用过以后，需要手动释放掉*w的内存空间529 * 调用示例：ret = gaussianWin(99,2.5, &w); 530 */531 dspErrorStatus gaussianWin(dspUint_16 N, dspDouble alpha, dspDouble **w)532 {533 dspUint_16 n;534 dspDouble k, beta, theta;535 dspDouble *ret;536 537 ret = (dspDouble *)malloc(N * sizeof(dspDouble));538 if(ret == NULL)539 return DSP_ERROR;540 541 for(n =0; n < N; n++)542 {543 if((N % 2) == 1)544 {545 k = n - (N - 1)/2;546 beta = 2 * alpha * (dspDouble)k / (N - 1); 547 }548 else549 {550 k = n - (N)/2;551 beta = 2 * alpha * (dspDouble)k / (N - 1); 552 }553 554 theta = pow(beta, 2);555 *(ret + n) = exp((-1) * (dspDouble)theta / 2);556 }557 558 *w = ret;559 560 return DSP_SUCESS;561 }562 563 /*564 *函数名：hammingWin565 *说明：计算hammingWin窗函数566 *输入：567 *输出：568 *返回：569 *调用：570 *其它：用过以后，需要手动释放掉*w的内存空间571 * 调用示例：ret = hammingWin(99, &w); 572 */573 dspErrorStatus hammingWin(dspUint_16 N, dspDouble **w)574 {575 dspUint_16 n;576 dspDouble *ret;577 ret = (dspDouble *)malloc(N * sizeof(dspDouble));578 if(ret == NULL)579 return DSP_ERROR;580 581 for(n = 0; n < N; n++)582 {583 *(ret + n) = 0.54 - 0.46 * cos (2 * PI * ( dspDouble )n / ( N - 1 ) );584 }585 586 *w = ret;587 588 return DSP_SUCESS;589 }590 591 /*592 *函数名：hannWin593 *说明：计算hannWin窗函数594 *输入：595 *输出：596 *返回：597 *调用：598 *其它：用过以后，需要手动释放掉*w的内存空间599 * 调用示例：ret = hannWin(99, &w); 600 */601 dspErrorStatus hannWin(dspUint_16 N, dspDouble **w)602 {603 dspUint_16 n;604 dspDouble *ret;605 ret = (dspDouble *)malloc(N * sizeof(dspDouble));606 if(ret == NULL)607 return DSP_ERROR;608 609 for(n = 0; n < N; n++)610 {611 *(ret + n) = 0.5 * ( 1 - cos( 2 * PI * (dspDouble)n / (N - 1)));612 }613 614 *w = ret;615 616 return DSP_SUCESS;617 }618 619 /*620 *函数名：kaiserWin621 *说明：计算kaiserWin窗函数622 *输入：623 *输出：624 *返回：625 *调用：besseli()第一类修正贝塞尔函数626 *其它：用过以后，需要手动释放掉*w的内存空间627 * 调用示例：ret = kaiserWin(99, 5, &w); 628 */629 dspErrorStatus kaiserWin(dspUint_16 N, dspDouble beta, dspDouble **w)630 {631 dspUint_16 n;632 dspDouble *ret;633 dspDouble theta;634 635 ret = (dspDouble *)malloc(N * sizeof(dspDouble));636 if(ret == NULL)637 return DSP_ERROR;638 639 for(n = 0; n < N; n++)640 {641 theta = beta * sqrt( 1 - pow( ( (2 * (dspDouble)n/(N -1)) - 1),2 ) );642 *(ret + n) = (dspDouble)besseli(0, theta, BESSELI_K_LENGTH) / besseli(0, beta, BESSELI_K_LENGTH);643 }644 645 *w = ret;646 647 return DSP_SUCESS;648 }649 650 /*651 *函数名：nuttalWin652 *说明：计算nuttalWin窗函数653 *输入：654 *输出：655 *返回：656 *调用：657 *其它：用过以后，需要手动释放掉*w的内存空间658 * 调用示例：ret = nuttalWin(99, &w); 659 */660 dspErrorStatus nuttalWin(dspUint_16 N, dspDouble **w)661 {662 dspUint_16 n;663 dspDouble *ret;664 665 ret = (dspDouble *)malloc(N * sizeof(dspDouble));666 if(ret == NULL)667 return DSP_ERROR;668 669 for(n = 0; n < N; n++)670 {671 *(ret + n) =NUTTALL_A0 - NUTTALL_A1 * cos(2 * PI * (dspDouble)n / (N - 1)) +\672 NUTTALL_A2 * cos(4 * PI * (dspDouble)n / (N - 1)) -\673 NUTTALL_A3 * cos(6 * PI * (dspDouble)n / (N - 1));674 675 }676 677 *w = ret;678 679 return DSP_SUCESS;680 }681 682 /*683 *函数名：parzenWin684 *说明：计算parzenWin窗函数685 *输入：686 *输出：687 *返回：688 *调用：689 *其它：用过以后，需要手动释放掉*w的内存空间690 * 调用示例：ret = parzenWin(99, &w); 691 */692 dspErrorStatus parzenWin(dspUint_16 N, dspDouble **w)693 {694 dspUint_16 n;695 dspDouble *ret;696 dspDouble alpha,k;697 698 ret = (dspDouble *)malloc(N * sizeof(dspDouble));699 if(ret == NULL)700 return DSP_ERROR;701 702 if(( N % 2) == 1)703 {704 for(n = 0; n < N; n++)705 {706 k = n - (N - 1) / 2;707 alpha = 2 * (dspDouble)myAbs(k) / N;708 if(myAbs(k) <= (N - 1) / 4)709 {710 *(ret + n) = 1 - 6 * pow(alpha,2) + 6 * pow(alpha, 3);711 }712 else713 {714 *(ret + n) = 2 * pow( (1 - alpha), 3 );715 }716 717 }718 }719 else720 {721 for(n = 0; n < N; n++)722 {723 k = n - (N - 1) / 2;724 alpha = 2 * (dspDouble)myAbs(k) / N;725 if(myAbs(k) <= (dspDouble)(N -1) / 4)726 {727 *(ret + n) = 1 - 6 * pow(alpha,2) + 6 * pow(alpha, 3);728 }729 else730 {731 *(ret + n) = 2 * pow( (1 - alpha), 3 );732 }733 734 }735 }736 737 738 739 *w = ret;740 741 return DSP_SUCESS;742 }743 744 /*745 *函数名：rectangularWin746 *说明：计算rectangularWin窗函数747 *输入：748 *输出：749 *返回：750 *调用：751 *其它：用过以后，需要手动释放掉*w的内存空间752 * 调用示例：ret = rectangularWin(99, &w); 753 */754 dspErrorStatus rectangularWin(dspUint_16 N, dspDouble **w)755 {756 dspUint_16 n;757 dspDouble *ret;758 759 ret = (dspDouble *)malloc(N * sizeof(dspDouble));760 if(ret == NULL)761 return DSP_ERROR;762 763 for(n = 0; n < N; n++)764 {765 *(ret + n) = 1; 766 }767 768 *w = ret;769 770 return DSP_SUCESS;771 }
欢迎多交流！
1 /*窗类型*/ 2 typedef enum 3 { 4 Bartlett = 0, 5 BartLettHann, 6 BlackMan, 7 BlackManHarris, 8 Bohman, 9 Chebyshev,10 FlatTop,11 Gaussian,12 Hamming,13 Hann,14 Kaiser,15 Nuttal,16 Parzen,17 Rectangular,18 Taylor, 19 Triangular,20 Tukey21 }winType;
别的不多说了，直接上干货。
1 /* 2 *file WindowFunction.h 3 *author Vincent Cui 4 *e-mail whcui1987@163.com 5 *version 0.3 6 *data 31-Oct-2014 7 *brief 各种窗函数的C语言实现 8 */ 9 10 11 #ifndef _WINDOWFUNCTION_H_12 #define _WINDOWFUNCTION_H_13 14 #include "GeneralConfig.h"15 16 #define BESSELI_K_LENGTH 1017 18 #define FLATTOPWIN_A0 0.21557899519 #define FLATTOPWIN_A1 0.4166315820 #define FLATTOPWIN_A2 0.27726315821 #define FLATTOPWIN_A3 0.08357894722 #define FLATTOPWIN_A4 0.00694736823 24 #define NUTTALL_A0 0.363581925 #define NUTTALL_A1 0.489177526 #define NUTTALL_A2 0.136599527 #define NUTTALL_A3 0.010641128 29 #define BLACKMANHARRIS_A0 0.3587530 #define BLACKMANHARRIS_A1 0.4882931 #define BLACKMANHARRIS_A2 0.1412832 #define BLACKMANHARRIS_A3 0.0116833 34 35 36 dspErrorStatus taylorWin(dspUint_16 N, dspUint_16 nbar, dspDouble sll, dspDouble **w);37 dspErrorStatus triangularWin(dspUint_16 N, dspDouble **w);38 dspErrorStatus tukeyWin(dspUint_16 N, dspDouble r, dspDouble **w);39 dspErrorStatus bartlettWin(dspUint_16 N, dspDouble **w);40 dspErrorStatus bartLettHannWin(dspUint_16 N, dspDouble **w);41 dspErrorStatus blackManWin(dspUint_16 N, dspDouble **w);42 dspErrorStatus blackManHarrisWin(dspUint_16 N, dspDouble **w);43 dspErrorStatus bohmanWin(dspUint_16 N, dspDouble **w);44 dspErrorStatus chebyshevWin(dspUint_16 N, dspDouble r, dspDouble **w);45 dspErrorStatus flatTopWin(dspUint_16 N, dspDouble **w);46 dspErrorStatus gaussianWin(dspUint_16 N, dspDouble alpha, dspDouble **w);47 dspErrorStatus hammingWin(dspUint_16 N, dspDouble **w);48 dspErrorStatus hannWin(dspUint_16 N, dspDouble **w);49 dspErrorStatus kaiserWin(dspUint_16 N, dspDouble beta, dspDouble **w);50 dspErrorStatus nuttalWin(dspUint_16 N, dspDouble **w);51 dspErrorStatus parzenWin(dspUint_16 N, dspDouble **w);52 dspErrorStatus rectangularWin(dspUint_16 N, dspDouble **w);53 54 55 56 #endif
1 /* 2 *file WindowFunction.c 3 *author Vincent Cui 4 *e-mail whcui1987@163.com 5 *version 0.3 6 *data 31-Oct-2014 7 *brief 各种窗函数的C语言实现 8 */ 9 10 11 #include "WindowFunction.h" 12 #include "GeneralConfig.h" 13 #include "MathReplenish.h" 14 #include "math.h" 15 #include <stdlib.h> 16 #include <stdio.h> 17 #include <string.h> 18 19 /*函数名：taylorWin 20 *说明：计算泰勒窗。泰勒加权函数 21 *输入： 22 *输出： 23 *返回： 24 *调用：prod()连乘函数 25 *其它：用过以后，需要手动释放掉*w的内存空间 26 * 调用示例：ret = taylorWin(99, 4, 40, &w); 注意此处的40是正数 表示-40dB 27 */ 28 dspErrorStatus taylorWin(dspUint_16 N, dspUint_16 nbar, dspDouble sll, dspDouble **w) 29 { 30 dspDouble A; 31 dspDouble *retDspDouble; 32 dspDouble *sf; 33 dspDouble *result; 34 dspDouble alpha,beta,theta; 35 dspUint_16 i,j; 36 37 /*A = R cosh(PI, A) = R*/ 38 A = (dspDouble)acosh(pow((dspDouble)10.0,(dspDouble)sll/20.0)) / PI; 39 A = A * A; 40 41 /*开出存放系数的空间*/ 42 retDspDouble = (dspDouble *)malloc(sizeof(dspDouble) * (nbar - 1)); 43 if(retDspDouble == NULL) 44 return DSP_ERROR; 45 sf = retDspDouble; 46 47 /*开出存放系数的空间*/ 48 retDspDouble = (dspDouble *)malloc(sizeof(dspDouble) * N); 49 if(retDspDouble == NULL) 50 return DSP_ERROR; 51 result = retDspDouble; 52 53 alpha = prod(1, 1, (nbar - 1)); 54 alpha *= alpha; 55 beta = (dspDouble)nbar / sqrt( A + pow((nbar - 0.5), 2) ); 56 for(i = 1; i <= (nbar - 1); i++) 57 { 58 *(sf + i - 1) = prod(1,1,(nbar -1 + i)) * prod(1,1,(nbar -1 - i)); 59 theta = 1; 60 for(j = 1; j <= (nbar - 1); j++) 61 { 62 theta *= 1 - (dspDouble)(i * i) / ( beta * beta * ( A + (j - 0.5) * (j - 0.5)) ); 63 } 64 *(sf + i - 1) = alpha * (dspDouble)theta / (*(sf + i - 1)); 65 } 66 67 /*奇数阶*/ 68 if((N % 2) == 1) 69 { 70 for(i = 0; i < N; i++) 71 { 72 alpha = 0; 73 for(j = 1; j <= (nbar - 1); j++) 74 { 75 alpha += (*(sf + j - 1)) * cos( 2 * PI * j * (dspDouble)(i - ((N-1)/2))/N ); 76 } 77 *(result + i) = 1 + 2 * alpha; 78 } 79 } 80 /*偶数阶*/ 81 else 82 { 83 for(i = 0; i < N; i++) 84 { 85 alpha = 0; 86 for(j = 1; j <= (nbar - 1); j++) 87 { 88 alpha += (*(sf + j - 1)) * cos( PI * j * (dspDouble)(2 * (i - (N/2)) + 1) / N ); 89 } 90 *(result + i) = 1 + 2 * alpha; 91 92 } 93 } 94 *w = result; 95 free(sf); 96 97 return DSP_SUCESS; 98 99 }100 101 102 /*103 *函数名：triangularWin104 *说明：计算三角窗函数105 *输入：106 *输出：107 *返回：108 *调用：109 *其它：用过以后，需要手动释放掉*w的内存空间110 * 调用示例：ret = triangularWin(99, &w); 111 */112 dspErrorStatus triangularWin(dspUint_16 N, dspDouble **w)113 {114 dspDouble *ptr;115 dspUint_16 i;116 117 ptr = (dspDouble *)malloc(N * sizeof(dspDouble));118 if(ptr == NULL)119 return DSP_ERROR;120 121 122 /*阶数为奇*/123 if((N % 2) == 1)124 {125 for(i = 0; i < ((N - 1)/2); i++)126 {127 *(ptr + i) = 2 * (dspDouble)(i + 1) / (N + 1);128 }129 for(i = ((N - 1)/2); i < N; i++)130 {131 *(ptr + i) = 2 * (dspDouble)(N - i) / (N + 1);132 }133 }134 /*阶数为偶*/135 else136 {137 for(i = 0; i < (N/2); i++)138 {139 *(ptr + i) = (i + i + 1) * (dspDouble)1 / N;140 }141 for(i = (N/2); i < N; i++)142 {143 *(ptr + i) = *(ptr + N - 1 - i);144 }145 }146 *w = ptr;147 148 return DSP_SUCESS;149 }150 151 /*152 *函数名：tukeyWin153 *说明：计算tukey窗函数154 *输入：155 *输出：156 *返回：linSpace()157 *调用：158 *其它：用过以后，需要手动释放掉*w的内存空间159 * 调用示例：ret = tukeyWin(99, 0.5, &w); 160 */161 dspErrorStatus tukeyWin(dspUint_16 N, dspDouble r, dspDouble **w)162 {163 dspErrorStatus retErrorStatus;164 dspUint_16 index;165 dspDouble *x,*result,*retPtr;166 dspDouble alpha;167 168 retErrorStatus = linSpace(0, 1, N, &x);169 if(retErrorStatus == DSP_ERROR)170 return DSP_ERROR;171 172 result = (dspDouble *)malloc(N * sizeof(dspDouble));173 if(result == NULL)174 return DSP_ERROR;175 176 /*r <= 0 就是矩形窗*/177 if(r <= 0)178 {179 retErrorStatus = rectangularWin(N, &retPtr);180 if(retErrorStatus == DSP_ERROR)181 return DSP_ERROR;182 /*将数据拷出来以后，释放调用的窗函数的空间*/183 memcpy(result, retPtr, ( N * sizeof(dspDouble)));184 free(retPtr);185 }186 /*r >= 1 就是汉宁窗*/187 else if(r >= 1)188 {189 retErrorStatus = hannWin(N, &retPtr);190 if(retErrorStatus == DSP_ERROR)191 return DSP_ERROR;192 /*将数据拷出来以后，释放调用的窗函数的空间*/193 memcpy(result, retPtr, ( N * sizeof(dspDouble)));194 free(retPtr);195 }196 else197 {198 for(index = 0; index < N; index++)199 {200 alpha = *(x + index);201 if(alpha < (r/2))202 {203 *(result + index) = (dspDouble)(1 + cos( 2 * PI * (dspDouble)(alpha - (dspDouble)r/2)/r))/2;204 }205 else if((alpha >= (r/2)) && (alpha <(1 - r/2)))206 {207 *(result + index) = 1;208 }209 else210 {211 *(result + index) = (dspDouble)(1 + cos( 2 * PI * (dspDouble)(alpha - 1 + (dspDouble)r/2)/r))/2;212 }213 214 }215 }216 217 free(x);218 219 *w = result;220 221 return DSP_SUCESS;222 223 }224 225 /*226 *函数名：bartlettWin227 *说明：计算bartlettWin窗函数228 *输入：229 *输出：230 *返回：231 *调用：232 *其它：用过以后，需要手动释放掉*w的内存空间233 * 调用示例：ret = bartlettWin(99, &w); 234 */235 dspErrorStatus bartlettWin(dspUint_16 N, dspDouble **w)236 {237 dspDouble *ret;238 dspUint_16 n;239 240 ret = (dspDouble *)malloc(N * sizeof(dspDouble));241 if(ret == NULL)242 return DSP_ERROR;243 244 for(n = 0; n < ( N - 1 ) / 2; n++)245 {246 *(ret + n) = 2 * (dspDouble)n / (N - 1);247 }248 249 for(n = ( N - 1 ) / 2; n < N; n++)250 {251 *(ret + n) = 2 - 2 * (dspDouble)n / (( N - 1 ));252 }253 254 *w = ret;255 256 return DSP_SUCESS;257 }258 259 /*260 *函数名：bartLettHannWin261 *说明：计算bartLettHannWin窗函数262 *输入：263 *输出：264 *返回：265 *调用：266 *其它：用过以后，需要手动释放掉*w的内存空间267 * 调用示例：ret = bartLettHannWin(99, &w); 268 */269 dspErrorStatus bartLettHannWin(dspUint_16 N, dspDouble **w)270 {271 dspUint_16 n;272 dspDouble *ret;273 274 ret = (dspDouble *)malloc(N * sizeof(dspDouble));275 if(ret == NULL)276 return DSP_ERROR;277 /*奇*/278 if(( N % 2 ) == 1)279 {280 for(n = 0; n < N; n++)281 {282 *(ret + n) = 0.62 - 0.48 * myAbs( ( (dspDouble)n / ( N - 1 ) ) - 0.5 ) + 0.38 * cos( 2 * PI * ( ((dspDouble)n / ( N - 1 ) ) - 0.5 ) );283 }284 for(n = 0; n < (N-1)/2; n++)285 {286 *(ret + n) = *(ret + N - 1 - n);287 }288 }289 /*偶*/290 else291 {292 for(n = 0; n < N; n++)293 {294 *(ret + n) = 0.62 - 0.48 * myAbs( ( (dspDouble)n / ( N - 1 ) ) - 0.5 ) + 0.38 * cos( 2 * PI * ( ((dspDouble)n / ( N - 1 ) ) - 0.5 ) );295 }296 for(n = 0; n < N/2; n++)297 {298 *(ret + n) = *(ret + N -1 - n);299 }300 }301 302 *w = ret;303 304 return DSP_SUCESS;305 306 }307 308 /*309 *函数名：blackManWin310 *说明：计算blackManWin窗函数311 *输入：312 *输出：313 *返回：314 *调用：315 *其它：用过以后，需要手动释放掉*w的内存空间316 * 调用示例：ret = blackManWin(99, &w); 317 */318 dspErrorStatus blackManWin(dspUint_16 N, dspDouble **w)319 {320 dspUint_16 n;321 dspDouble *ret;322 ret = (dspDouble *)malloc(N * sizeof(dspDouble));323 if(ret == NULL)324 return DSP_ERROR;325 326 for(n = 0; n < N; n++)327 {328 *(ret + n) = 0.42 - 0.5 * cos(2 * PI * (dspDouble)n / ( N - 1 )) + 0.08 * cos( 4 * PI * ( dspDouble )n / ( N - 1 ) );329 }330 331 *w = ret;332 333 return DSP_SUCESS;334 }335 336 /*337 *函数名：blackManHarrisWin338 *说明：计算blackManHarrisWin窗函数339 *输入：340 *输出：341 *返回：342 *调用：343 *其它：用过以后，需要手动释放掉*w的内存空间344 * 调用示例：ret = blackManHarrisWin(99, &w); 345 * minimum 4-term Blackman-harris window -- From Matlab346 */347 dspErrorStatus blackManHarrisWin(dspUint_16 N, dspDouble **w)348 {349 dspUint_16 n;350 dspDouble *ret;351 352 ret = (dspDouble *)malloc(N * sizeof(dspDouble));353 if(ret == NULL)354 return DSP_ERROR;355 356 for(n = 0; n < N; n++)357 {358 *(ret + n) = BLACKMANHARRIS_A0 - BLACKMANHARRIS_A1 * cos( 2 * PI * (dspDouble)n / (N) ) + \359 BLACKMANHARRIS_A2 * cos( 4 * PI * (dspDouble)n/ (N) ) - \360 BLACKMANHARRIS_A3 * cos( 6 * PI * (dspDouble)n/ (N) );361 }362 363 *w = ret;364 365 return DSP_SUCESS;366 }367 368 /*369 *函数名：bohmanWin370 *说明：计算bohmanWin窗函数371 *输入：372 *输出：373 *返回：374 *调用：375 *其它：用过以后，需要手动释放掉*w的内存空间376 * 调用示例：ret = bohmanWin(99, &w); 377 */378 dspErrorStatus bohmanWin(dspUint_16 N, dspDouble **w)379 {380 dspUint_16 n;381 dspDouble *ret;382 dspDouble x;383 ret = (dspDouble *)malloc(N * sizeof(dspDouble));384 if(ret == NULL)385 return DSP_ERROR;386 387 for(n = 0; n < N; n++)388 {389 x = -1 + n * (dspDouble)2 / ( N - 1 ) ;390 /*取绝对值*/391 x = x >= 0 ? x : ( x * ( -1 ) );392 *(ret + n) = ( 1 - x ) * cos( PI * x) + (dspDouble)(1 / PI) * sin( PI * x);393 }394 395 *w = ret;396 397 return DSP_SUCESS;398 }399 400 /*401 *函数名：chebyshevWin402 *说明：计算chebyshevWin窗函数403 *输入：404 *输出：405 *返回：406 *调用：407 *其它：用过以后，需要手动释放掉*w的内存空间408 * 调用示例：ret = chebyshevWin(99,100, &w); 409 */410 dspErrorStatus chebyshevWin(dspUint_16 N, dspDouble r, dspDouble **w)411 {412 dspUint_16 n,index;413 dspDouble *ret;414 dspDouble x, alpha, beta, theta, gama;415 416 ret = (dspDouble *)malloc(N * sizeof(dspDouble));417 if(ret == NULL)418 return DSP_ERROR;419 420 421 /*10^(r/20)*/422 theta = pow((dspDouble)10, (dspDouble)(myAbs(r)/20));423 beta = pow(cosh(acosh(theta)/(N - 1)),2);424 alpha = 1 - (dspDouble)1 / beta;425 426 if((N % 2) == 1)427 {428 /*计算一半的区间*/429 for( n = 1; n < ( N + 1 ) / 2; n++ )430 {431 gama = 1;432 for(index = 1; index < n; index++)433 {434 x = index * (dspDouble)( N - 1 - 2 * n + index) /(( n - index ) * (n + 1 -index));435 gama = gama * alpha * x + 1;436 }437 *(ret + n) = (N - 1) * alpha * gama; 438 }439 440 theta = *( ret + (N - 1)/2 );441 *ret = 1;442 443 for(n = 0; n < ( N + 1 ) / 2; n++ )444 {445 *(ret + n) = (dspDouble)(*(ret + n)) / theta;446 }447 448 /*填充另一半*/449 for(; n < N; n++)450 {451 *(ret + n) = ret[N - n - 1];452 }453 }454 else455 {456 /*计算一半的区间*/457 for( n = 1; n < ( N + 1 ) / 2; n++ )458 {459 gama = 1;460 for(index = 1; index < n; index++)461 {462 x = index * (dspDouble)( N - 1 - 2 * n + index) /(( n - index ) * (n + 1 -index));463 gama = gama * alpha * x + 1;464 }465 *(ret + n) = (N - 1) * alpha * gama; 466 }467 468 theta = *( ret + (N/2) - 1);469 *ret = 1;470 471 for(n = 0; n < ( N + 1 ) / 2; n++ )472 {473 *(ret + n) = (dspDouble)(*(ret + n)) / theta;474 }475 476 /*填充另一半*/477 for(; n < N; n++)478 {479 *(ret + n) = ret[N - n - 1];480 }481 }482 483 484 *w = ret;485 486 return DSP_SUCESS;487 }488 489 /*490 *函数名：flatTopWin491 *说明：计算flatTopWin窗函数492 *输入：493 *输出：494 *返回：495 *调用：496 *其它：用过以后，需要手动释放掉*w的内存空间497 * 调用示例：ret = flatTopWin(99, &w); 498 */499 dspErrorStatus flatTopWin(dspUint_16 N, dspDouble **w)500 {501 dspUint_16 n;502 dspDouble *ret;503 ret = (dspDouble *)malloc(N * sizeof(dspDouble));504 if(ret == NULL)505 return DSP_ERROR;506 507 for(n = 0; n < N; n++)508 {509 *(ret + n) = FLATTOPWIN_A0 - FLATTOPWIN_A1 * cos(2 * PI * (dspDouble)n / (N - 1)) +\510 FLATTOPWIN_A2 * cos(4 * PI * (dspDouble)n / (N - 1)) -\511 FLATTOPWIN_A3 * cos(6 * PI * (dspDouble)n / (N - 1)) +\512 FLATTOPWIN_A4 * cos(8 * PI * (dspDouble)n / (N - 1));513 }514 515 *w = ret;516 517 return DSP_SUCESS;518 }519 520 521 /*522 *函数名：gaussianWin523 *说明：计算gaussianWin窗函数524 *输入：525 *输出：526 *返回：527 *调用：528 *其它：用过以后，需要手动释放掉*w的内存空间529 * 调用示例：ret = gaussianWin(99,2.5, &w); 530 */531 dspErrorStatus gaussianWin(dspUint_16 N, dspDouble alpha, dspDouble **w)532 {533 dspUint_16 n;534 dspDouble k, beta, theta;535 dspDouble *ret;536 537 ret = (dspDouble *)malloc(N * sizeof(dspDouble));538 if(ret == NULL)539 return DSP_ERROR;540 541 for(n =0; n < N; n++)542 {543 if((N % 2) == 1)544 {545 k = n - (N - 1)/2;546 beta = 2 * alpha * (dspDouble)k / (N - 1); 547 }548 else549 {550 k = n - (N)/2;551 beta = 2 * alpha * (dspDouble)k / (N - 1); 552 }553 554 theta = pow(beta, 2);555 *(ret + n) = exp((-1) * (dspDouble)theta / 2);556 }557 558 *w = ret;559 560 return DSP_SUCESS;561 }562 563 /*564 *函数名：hammingWin565 *说明：计算hammingWin窗函数566 *输入：567 *输出：568 *返回：569 *调用：570 *其它：用过以后，需要手动释放掉*w的内存空间571 * 调用示例：ret = hammingWin(99, &w); 572 */573 dspErrorStatus hammingWin(dspUint_16 N, dspDouble **w)574 {575 dspUint_16 n;576 dspDouble *ret;577 ret = (dspDouble *)malloc(N * sizeof(dspDouble));578 if(ret == NULL)579 return DSP_ERROR;580 581 for(n = 0; n < N; n++)582 {583 *(ret + n) = 0.54 - 0.46 * cos (2 * PI * ( dspDouble )n / ( N - 1 ) );584 }585 586 *w = ret;587 588 return DSP_SUCESS;589 }590 591 /*592 *函数名：hannWin593 *说明：计算hannWin窗函数594 *输入：595 *输出：596 *返回：597 *调用：598 *其它：用过以后，需要手动释放掉*w的内存空间599 * 调用示例：ret = hannWin(99, &w); 600 */601 dspErrorStatus hannWin(dspUint_16 N, dspDouble **w)602 {603 dspUint_16 n;604 dspDouble *ret;605 ret = (dspDouble *)malloc(N * sizeof(dspDouble));606 if(ret == NULL)607 return DSP_ERROR;608 609 for(n = 0; n < N; n++)610 {611 *(ret + n) = 0.5 * ( 1 - cos( 2 * PI * (dspDouble)n / (N - 1)));612 }613 614 *w = ret;615 616 return DSP_SUCESS;617 }618 619 /*620 *函数名：kaiserWin621 *说明：计算kaiserWin窗函数622 *输入：623 *输出：624 *返回：625 *调用：besseli()第一类修正贝塞尔函数626 *其它：用过以后，需要手动释放掉*w的内存空间627 * 调用示例：ret = kaiserWin(99, 5, &w); 628 */629 dspErrorStatus kaiserWin(dspUint_16 N, dspDouble beta, dspDouble **w)630 {631 dspUint_16 n;632 dspDouble *ret;633 dspDouble theta;634 635 ret = (dspDouble *)malloc(N * sizeof(dspDouble));636 if(ret == NULL)637 return DSP_ERROR;638 639 for(n = 0; n < N; n++)640 {641 theta = beta * sqrt( 1 - pow( ( (2 * (dspDouble)n/(N -1)) - 1),2 ) );642 *(ret + n) = (dspDouble)besseli(0, theta, BESSELI_K_LENGTH) / besseli(0, beta, BESSELI_K_LENGTH);643 }644 645 *w = ret;646 647 return DSP_SUCESS;648 }649 650 /*651 *函数名：nuttalWin652 *说明：计算nuttalWin窗函数653 *输入：654 *输出：655 *返回：656 *调用：657 *其它：用过以后，需要手动释放掉*w的内存空间658 * 调用示例：ret = nuttalWin(99, &w); 659 */660 dspErrorStatus nuttalWin(dspUint_16 N, dspDouble **w)661 {662 dspUint_16 n;663 dspDouble *ret;664 665 ret = (dspDouble *)malloc(N * sizeof(dspDouble));666 if(ret == NULL)667 return DSP_ERROR;668 669 for(n = 0; n < N; n++)670 {671 *(ret + n) =NUTTALL_A0 - NUTTALL_A1 * cos(2 * PI * (dspDouble)n / (N - 1)) +\672 NUTTALL_A2 * cos(4 * PI * (dspDouble)n / (N - 1)) -\673 NUTTALL_A3 * cos(6 * PI * (dspDouble)n / (N - 1));674 675 }676 677 *w = ret;678 679 return DSP_SUCESS;680 }681 682 /*683 *函数名：parzenWin684 *说明：计算parzenWin窗函数685 *输入：686 *输出：687 *返回：688 *调用：689 *其它：用过以后，需要手动释放掉*w的内存空间690 * 调用示例：ret = parzenWin(99, &w); 691 */692 dspErrorStatus parzenWin(dspUint_16 N, dspDouble **w)693 {694 dspUint_16 n;695 dspDouble *ret;696 dspDouble alpha,k;697 698 ret = (dspDouble *)malloc(N * sizeof(dspDouble));699 if(ret == NULL)700 return DSP_ERROR;701 702 if(( N % 2) == 1)703 {704 for(n = 0; n < N; n++)705 {706 k = n - (N - 1) / 2;707 alpha = 2 * (dspDouble)myAbs(k) / N;708 if(myAbs(k) <= (N - 1) / 4)709 {710 *(ret + n) = 1 - 6 * pow(alpha,2) + 6 * pow(alpha, 3);711 }712 else713 {714 *(ret + n) = 2 * pow( (1 - alpha), 3 );715 }716 717 }718 }719 else720 {721 for(n = 0; n < N; n++)722 {723 k = n - (N - 1) / 2;724 alpha = 2 * (dspDouble)myAbs(k) / N;725 if(myAbs(k) <= (dspDouble)(N -1) / 4)726 {727 *(ret + n) = 1 - 6 * pow(alpha,2) + 6 * pow(alpha, 3);728 }729 else730 {731 *(ret + n) = 2 * pow( (1 - alpha), 3 );732 }733 734 }735 }736 737 738 739 *w = ret;740 741 return DSP_SUCESS;742 }743 744 /*745 *函数名：rectangularWin746 *说明：计算rectangularWin窗函数747 *输入：748 *输出：749 *返回：750 *调用：751 *其它：用过以后，需要手动释放掉*w的内存空间752 * 调用示例：ret = rectangularWin(99, &w); 753 */754 dspErrorStatus rectangularWin(dspUint_16 N, dspDouble **w)755 {756 dspUint_16 n;757 dspDouble *ret;758 759 ret = (dspDouble *)malloc(N * sizeof(dspDouble));760 if(ret == NULL)761 return DSP_ERROR;762 763 for(n = 0; n < N; n++)764 {765 *(ret + n) = 1; 766 }767 768 *w = ret;769 770 return DSP_SUCESS;771 }
欢迎多交流！
阅读全文
0 0
- 窗函数的C语言实现
- 窗函数的C语言实现
- C语言可变参数函数的实现
- C语言可变参数函数的实现
- C语言实现的画直方图函数
- trim函数的C语言实现
- C语言可变参数的函数实现
- C语言实现：函数的递归调用
- C语言atoi函数的实现
- C语言实现的trim函数str_trim
- 内存拷贝函数的C语言实现
- c语言字符串函数的实现
- C语言函数调用的汇编实现
- c语言 分割函数的实现
- trim函数的C语言实现
- C语言是实现函数的重载
- C语言可变参数函数的实现
- 【C语言】strcpy函数的实现
- 几种数据库的大数据批量插入
- SpringBoot报错（六）热部署报错
- 云原生落地蓝图：Educate、Enable、Enforce
- flume配置启动
- tensorflow里面共享变量、name_scope, variable_scope等如何理解
- 窗函数的C语言实现
- 概率算法中的Monte carlo算法
- h5开发之于APP有哪些技术上的革新？
- js中的闭包
- 2017/10/16 弹层组件文档
- 新人学习笔记-分布式基础学习-分布式文件系统
- jstorm、storm 常见问题
- 网狐棋牌 权位过滤
- TCP连接
